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1.  Introduction 


The  split  Hopkinson  pressure  bar  (SHPB),  also  known  as  the  Kolsky  bar,  is  an  experimental 
technique  that  is  widely  used  for  studying  the  strain  rate  sensitivity  of  inelastic  materials  in  a 
state  of  compressive  uniaxial  stress.  The  standard  relations  for  analyzing  data  from  SHPB  tests 
rely  on  the  assumption  of  uniform  stress,  strain  and  strain  rate  throughout  the  specimen.  These 
conditions,  as  well  as  a  nearly  constant  nominal  strain  rate,  can  often  be  achieved  after  an  initial 
“ringing-up”  period.  The  review  articles  by  Gray  (1)  and  Gama  et  al.  (2)  are  recommended  for 
historical  background,  discussions  of  experimental  procedures  and  data  analysis  techniques,  and 
additional  references  to  the  literature. 

The  most  common  version  of  the  SHPB  test  assumes  that  the  two  pressure  bars  undergo  small 
strain,  linear  elastic  deformations,  that  the  specimen  remains  in  contact  with  the  pressure  bars, 
and  that  the  faces  of  the  bars  in  contact  with  the  specimen  remain  planar.  These  conditions 
require  that  the  specimen  be  softer  than  the  pressure  bars.  However,  difficulties  arise  in  the 
analysis  of  the  experimental  data  when  the  specimen  is  extremely  soft  relative  to  the  bars.  These 
difficulties  as  well  as  techniques  that  have  been  developed  to  overcome  them  are  discussed  in 
Gray  ( 1 )  and  Gama  et  al.  (2)  and  also  in  Gray  and  Blumenthal  (i),  Song  and  Chen  ( 4 ),  Moy  et  al. 
(5),  and  Song  et  al.  (<5).  One  source  of  difficulty  is  that  the  signal  in  the  transmission  bar  may  be 
too  weak  to  provide  accurate  measurements  of  stress  in  the  specimen;  this  issue  is  not  addressed 
here.  Another  difficulty  when  the  imposed  axial  strain  rate  is  sufficiently  high  is  that  the  stress 
and  strain  within  a  soft  specimen  may  be  non-uniform  and  the  stress  state  non-uniaxial 
throughout  the  test.  Such  conditions  either  invalidate  the  test  or  at  least  require  appropriate 
“inertial”  corrections  to  the  test  data. 

The  original  objective  of  the  present  work  was  to  perform  numerical  simulations  on  very  soft 
materials  in  order  to  test  the  validity  of  some  approximate,  analytical,  inertial  corrections  for 
SHPB  test  data.  These  relations,  which  were  developed  by  Scheidler  (7),  utilized  the  simplifying 
approximation  of  constant  volume  and  thus  are  limited  to  nearly  incompressible  materials,  that 
is,  to  materials  that  are  substantially  stiffer  in  dilatation  than  in  shear.  Examples  of  soft,  nearly 
incompressible  materials  include  many  types  of  rubber  as  well  as  many  biological  materials  due 
to  their  high  water  content.  Ballistic  gelatin,  which  is  80-90%  water  by  mass,  also  falls  into  this 
category;  it  is  widely  used  as  a  tissue  surrogate  in  impact  and  penetration  tests. 

We  performed  numerous  axisymmetric  finite-element  simulations  of  SHPB  tests  on  soft,  nearly 
incompressible  materials  using  the  commercial  code  LS-DYNA  (8).  The  initial  simulations  used 
a  nonlinear  elastic  model  for  the  specimen,  namely,  LS-DYNA’s  compressible  version  of  the 
Mooney-Rivlin  model  (8,  9).  The  model  was  calibrated  to  give  rough  agreement  with  the 
uniaxial  compression  data  on  ballistic  gelatin  obtained  by  Moy  et  al.  (5)  and  with  the  pressure- 
volume  data  on  ballistic  gelatin  in  Aihaiti  and  Hemley  (10).  We  used  a  relatively  thin  specimen 
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(1.45-mm  thickness,  12.7-mm  diameter)  and  a  nominal  strain  rate  of  2500/s  in  all  simulations. 
These  were  chosen  to  agree  with  the  test  conditions  in  Moy  et  al.  (5).  For  these  initial 
simulations  we  used  solid  (as  opposed  to  annular)  specimens. 

Much  to  our  surprise,  we  observed  that  small  gaps  formed  along  both  the  specimen-incident  bar 
(S-IB)  interface  and  the  specimen-transmission  bar  (S-TB)  interface.  The  size  of  these  gaps 
ranged  from  sub-micron  to  as  large  as  43  pm,  depending  on  the  rise  time  of  the  loading  wave  in 
the  incident  bar,  the  time  after  the  arrival  of  the  wave  at  the  specimen,  the  particular  interface, 
and  the  radial  location  on  that  interface.  The  gaps  closed  and  re-opened  multiple  times,  and  did 
not  necessarily  form  over  the  entire  face  of  the  specimen.  However,  for  the  S-IB  interface  and 
the  shortest  rise  time,  a  gap  existed  over  most  of  the  face  of  the  specimen  for  axial  strains 
ranging  from  4%  to  19%,  and  out  to  half  of  the  specimen  radius  for  strains  up  to  28%. 1  The 
faces  of  the  pressure  bars  remained  very  nearly  planar,  as  expected;  the  gaps  formed  as  the 
specimen  moved  axially  inward  toward  its  center.  Large  pressure  spikes  were  observed  on  the 
axis  of  the  bars  and  the  specimen  as  these  gaps  closed  up  and  the  specimen  “slapped”  the  bars. 

To  the  best  of  our  knowledge,  this  phenomenon  has  not  been  reported  in  either  the  experimental 
or  the  computational  literature  on  SHPB  tests,2  although  certain  features  of  previously  reported 
experimental  data  appear  to  be  consistent  with  the  opening  and  closing  of  gaps. 

The  assumption  of  contact  between  the  specimen  and  the  bars  is  implicit  in  any  analysis  of  the 
SHPB  test,  including  the  inertial  corrections  that  we  had  intended  to  validate.  Consequently,  the 
focus  of  our  study  shifted  to  the  gap  phenomenon.  In  particular,  we  sought  to  determine  whether 
the  gaps  were  a  numerical  artifact.  The  dependence  of  gap  formation  on  the  following 
conditions  was  examined: 

1.  Constitutive  Model  for  the  Specimen: 

a)  Model  type: 

i)  nonlinear  elastic  (compressible  Mooney-Rivlin) 

ii)  linear  elastic 

b)  Model  parameters  (linear  elastic  case  only) 

2.  Loading  Condition: 

a)  Linear  ramp  (1  or  25  ps)  to  constant  velocity  at  the  far  end  of  the  incident  bar 

b)  Impact  of  a  striker  bar  on  the  far  end  of  a  long  incident  bar 

c)  Direct  impact  of  a  striker  bar  on  the  specimen 


1 

1  Gaps  at  the  S-IB  interface  are  clearly  visible  in  figure  13  (solid  specimen)  and  figure  C-l  (annular  specimen). 

2  However,  a  summary  of  some  of  the  results  in  this  report  is  given  the  proceedings  paper  (11). 
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3.  Geometry: 

a)  Pressure  bar  diameter 

b)  Specimen  geometry: 

i)  solid 

ii)  annular 

4.  Boundary  Condition  on  the  Lateral  Surfaces: 

a)  Stress-free 

b)  Radially  constrained 

5.  Computational  Parameters: 

a)  Specimen  mesh  size 

b)  Contact  algorithm  parameters 

6.  The  code  used  for  the  simulations 

The  organization  of  the  report  can  be  inferred  from  the  table  of  contents.  Sections  2,  3,  4,  5.1, 
5.2,  and  appendix  A  contain  preliminary  and  background  material.  The  results  of  a  selected  set 
of  numerical  simulations  are  presented  in  detail  in  sections  5.3,  6,  7,  8,  and  appendix  B. 
Appendix  C  consists  of  less  detailed  summaries  of  simulations  for  which  some  of  the  conditions 
differed  from  those  in  the  main  body  of  the  report.  These  conditions  are  italicized  in  the  list 
above.  In  particular,  all  of  the  numerical  simulations  discussed  in  the  main  body  of  this  report 
were  performed  with  LS-DYNA.  In  appendix  C-3,  we  summarize  the  results  of  some  analogous 
simulations  performed  by  Bryan  Love  at  the  U.S.  Army  Research  Laboratory  (ARL)  using  the 
finite  element  code  PRESTO  (from  Sandia  National  Laboratories)  with  loading  condition 
described  in  2.b  above. 


2.  Problem  Description 


2.1  Problem  Geometry 

The  SHPB  experiment  was  modeled  using  the  three  parts  shown  in  figure  1,  which  is  not  to 
scale.  As  indicated  in  the  figure,  a  cylindrical  coordinate  system  was  used  with  the  origin 
located  on  the  axis  of  symmetry  or  centerline  of  the  specimen  and  the  pressure  bars  at  the  S-IB 
interface.  The  reference  or  material  coordinates  ( R ,  Z)  are  the  coordinates  in  the  undeformed 
state  at  time  t  =  0;  the  corresponding  coordinates  in  the  deformed  state  are  denoted  by  (r,  z).  For 
axisymmetric  deformations,  the  Z-axis,  the  z-axis,  the  set  of  points  with  r  =  0,  and  the  set  of 
points  with  R  =  0  all  describe  the  centerline.  The  positive  direction  of  the  Z-axis  is  along  the 
direction  of  propagation  of  the  initial  loading  wave,  that  is,  from  incident  bar  to  specimen  to 
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transmission  bar.  All  coordinate  values  and  problem  dimensions  are  given  in  millimeters, 
whereas  element  sizes  (section  4.1)  and  gap  sizes  are  specified  in  microns. 
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Figure  1.  Geometry  and  coordinate  system  for  the  SHPB  simulations  (not  drawn  to  scale). 


A  solid,  disk-shaped  specimen  (S)  was  used  for  most  of  the  simulations  described  in  this  report. 
A  brief  summary  of  simulations  involving  annular,  washer-shaped  specimens  is  given  in 
appendices  C-2  and  C-3.  Unless  specified  otherwise,  a  solid  specimen  is  assumed  in  subsequent 
discussions.  The  specimen  has  initial  thickness  (or  length)  Ls  =  1.45  mm  and  an  initial  diameter 
of  12.7  mm,  giving  a  length-to-diameter  ratio  of  0.1 14.  Relatively  thin  specimens  are  commonly 
used  to  decrease  the  time  it  takes  to  “ring-up”  to  a  uniform  state. 

The  incident  bar  (IB)  and  the  transmission  bar  (TB)  have  the  same  diameter.  Two  bar  diameters 
were  considered,  25.6  mm  and  19.0  mm.  Most  of  the  simulations  were  done  with  the  larger  bar 
diameter.  Discussion  of  simulations  with  the  smaller  bar  diameter  is  confined  to  section  6.4  and 
appendix  C-2.  The  length  of  the  incident  bar  used  in  the  simulations  is  Lib  =  768  mm.  This 
gives  a  length-to-diameter  ratio  of  30.0  for  the  larger  diameter  bar  and  40.4  for  the  smaller 
diameter  bar.  The  length  of  the  transmission  bar  used  in  the  simulations  is  LTB  =  256  mm,  which 
is  1/3  the  length  of  the  incident  bar  and  gives  length-to-diameter  ratios  of  10.0  and  13.5. 

The  dimensions  of  the  specimen  and  the  smaller  diameter  bar  were  chosen  to  agree  with  those 
used  in  the  experimental  study  on  gelatin  by  Moy  et  al.  (5).  However,  we  used  the  larger 
diameter  bar  in  most  of  our  simulations.  This  was  motivated  by  our  original  goal  of  studying 
inertial  effects  at  large  strains.  We  wished  to  obtain  nominal  axial  strains  of  up  to  75%  in 
compression,  which  corresponds  to  an  axial  stretch  of  0.25,  that  is,  a  specimen  compressed  to  !4 
of  its  original  thickness.  This  results  in  a  radial  stretch  of  2.0,  assuming  incompressibility  and  a 
homogeneous  deformation.  Hence  the  diameter  of  the  specimen  would  increase  by  a  factor  of  2, 
that  is,  from  12.7  to  25.4  mm.  And  because  the  deformed  specimen  must  not  extend  beyond  the 
pressure  bars,  we  chose  a  diameter  of  25.6  mm  for  the  bars. 
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We  will  often  refer  to  the  specimen  and  bar  radii  rather  than  the  diameters.  The  initial  radius  of 
the  specimen  is  Rs  =  6.35  mm,  and  the  two  pressure  bar  radii  considered  are  Rb  =  12.8  mm  and 
Rb  =  9.5  mm. 


2.2  Loading  and  Boundary  Conditions 

In  the  simplest  version  of  the  SHPB  test,  a  striker  bar  impacts  the  incident  bar  and  generates  a 
compressive  stress  pulse  which  propagates  down  the  length  of  the  incident  bar  and  subsequently 
into  the  specimen  and  the  transmission  bar.  A  direct  impact  of  a  striker  bar  on  the  incident  bar 
would  have  been  easy  to  simulate,  but  it  can  also  be  approximated  reasonably  well  by  a 
prescribed  step  in  the  axial  velocity  at  the  end  of  the  incident  bar. 


For  the  purpose  of  reducing  inertial  effects  in  soft  specimens,  it  is  beneficial  to  insert  a  pulse 
shaper  between  the  striker  and  incident  bars  (4—6).  This  is  a  thin  disk  of  another  soft  material 
(generally  not  the  same  as  the  specimen)  which  serves  to  smooth  the  incident  pulse  and  increase 
its  rise  time.  To  avoid  the  complexity  of  modeling  the  pulse  shaping  material  and,  in  particular, 
the  computational  issues  associated  with  the  extremely  large  deformations  of  the  pulse  shaper, 
we  chose  instead  to  approximate  the  effect  the  pulse  shaper  on  the  incident  pulse.  As  in  the 
computational  study  by  Song  et  al.  (6),  we  prescribed  an  axial  velocity,  v-(j),  on  the  impact  end 
of  the  incident  bar.  This  prescribed  velocity  has  an  initial  value  of  zero  and  rises  linearly  with 
time  until  attaining  a  plateau  value  v0  at  time  tR\ 


v,(0  = 


v—  ;  o <t<tR 

Vn  ;  t>tB 


(1) 


We  set  Vo  to  1.8125  m/s  in  all  calculations.  This  value  was  chosen  to  produce  a  nominal  axial 
strain  rate  of  approximately  2500/s  in  the  specimen,  in  agreement  with  the  tests  in  Moy  et  al.  (5). 
We  considered  two  values  for  the  initial  rise  time  tR,  1  and  25  ps. 

In  an  SHPB  test,  the  lateral  surfaces  of  the  pressure  bars  (R  =  Rb)  and  the  specimen  (R  =  Rs)  are 
unconstrained  and  hence  stress-free.  The  same  condition  was  used  in  the  numerical  simulations 
with  two  exceptions.  For  these  two  exceptional  cases  (discussed  in  appendix  C-l.l),  the  lateral 
surfaces  were  constrained  in  order  to  study  the  effects  of  eliminating  the  unloading  waves  from 
these  surfaces  and  reducing  the  radial  acceleration  in  the  specimen. 


Friction  between  the  faces  of  the  specimen  and  the  bars  can  result  in  bulging  of  the  specimen  on 
compression.  In  SHPB  tests  the  specimen  faces  are  lubricated  to  eliminate  (or  at  least 
substantially  reduce)  the  friction.  Consequently,  in  the  numerical  simulations  the  S-IB  and  S-TB 
interfaces  were  modeled  with  frictionless  contact  (see  section  4.2). 


At  the  far  end  of  the  transmission  bar,  a  non-reflecting  boundary  condition  was  used;  refer  to 
sections  4.2  and  5.2  for  further  discussion. 
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2.3  Material  Properties 


For  soft  specimens,  pressure  bar  materials  with  low  mechanical  impedance  are  often  used  to 
improve  the  signal  to  noise  ratio  (1-6).  For  example,  the  tests  reported  in  Moy  et  al.  (5)  and 
Song  et  al.  (6)  utilized  7075-T6  aluminum  bars.  For  the  imposed  velocity  history  discussed 
above,  the  stress  in  an  aluminum  bar  remains  well  below  yield  (see  section  5.1).  We  used  an 
isotropic,  linear  elastic  constitutive  model  for  the  incident  and  transmission  bars,  with  material 
properties  appropriate  for  aluminum.  These  properties  (obtained  from  Song  et  al.  [(5])  are  a 
density  of  2.7  g/cm  ,  a  Young’s  modulus  of  68  GPa,  and  a  Poisson’s  ratio  of  0.33. 

As  mentioned  in  the  Introduction,  we  chose  ballistic  gelatin  as  a  representative  soft,  nearly 
incompressible  material,  and  for  this  we  used  a  compressible  version  of  the  classical  Mooney- 
Rivlin  constitutive  model  for  incompressible  materials.  This  isotropic,  nonlinear  elastic 
constitutive  model  and  its  calibration  are  discussed  briefly  in  section  3  and  in  more  detail  in 
appendix  A,  and  the  results  of  simulations  with  this  model  are  given  in  sections  6  and  7  and 
appendices  B  and  C.  At  this  point  we  simply  note  that  the  density  of  the  specimen  was  chosen  to 
be  that  of  water,  1  g/cm' .  This  same  density  was  also  used  in  subsequent  simulations  with  a 
linear  elastic  model  for  the  specimen  (see  section  8). 3  These  simulations  were  performed  to 
determine  whether  or  not  the  formation  of  gaps  was  attributable  to  the  nonlinearity  and/or  near 
incompressibility  of  the  constitutive  model.  It  turned  out  that  gaps  formed  for  all  but  one  of  the 
six  sets  of  material  parameters  used  in  the  linear  elastic  model  for  the  specimen. 

Table  1  lists  the  density  and  linear  elastic  properties  used  for  the  pressure  bars  and  the  specimen 
in  the  various  simulations.  For  the  ballistic  gelatin  specimen,  these  parameters  represent  the 
small  strain,  linear  elastic  approximation  to  the  compressible  Mooney-Rivlin  model  (section  3) 
and  so  do  not  accurately  reflect  the  large  strain  response.  Nevertheless,  the  table  allows  one  to 
quickly  assess  the  relative  stiffnesses  of  the  materials.  E,  G,  K,  and  L  denote  the  Young’s 
modulus,  shear  modulus,  bulk  modulus,  and  longitudinal  modulus;  they  govern  the  response  to 
uniaxial  stress,  simple  shear,  pure  dilatation,  and  uniaxial  strain,  respectively,  v  and  p  denote 
Poisson’s  ratio  and  the  initial  density.  Since  all  the  materials  are  isotropic,  any  two  of  the 
material  parameters  E,  G,  K,  L,  and  v  determine  the  other  three.  For  example, 

2  1+v  l—v  4 

K= - G  and  L  = - - - E  =  K  +  -G.  (2) 

3  l-2v  (l  +  v)(l-2v)  3 


3  However,  the  linear  elastic  constants  used  in  these  simulations  are  not  necessarily  those  appropriate  for  the  small  strain 
response  of  ballistic  gelatin. 
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Table  1.  Linear  elastic  constants  for  materials  used  in  the  numerical  simulations. 


Material 

Density 

(g/cm3) 

Elastic  Moduli 
(GPa) 

Dimensionless  Elastic 
Constants 

Wave  Speeds 
(nini/ps) 

Impedance  Ratios 
(specimen/bar) 

P 

E 

G 

K 

L 

V 

G/K 

cL 

Ce 

CG 

Based 

on  cL 

Based 

on  cE 

Aluminum 

Bars 

2.70 

68.0 

25.6 

66.7 

101. 

0.33 

0.383 

6.11 

5.02 

3.08 

X 

X 

Ballistic 

Gelatin 

Specimen 

1.00 

2.40  xlCT4 

8.00  xl0~5 

4.00 

4.00 

0.49999 

2.00  xlCT5 

2.00 

0.0155 

0.00894 

0.121 

0.00114 

7.17 

7.17 

0.49999 

2.00  xl0~5 

2.68 

0.162 

1.43x10"* 

0.717 

0.717 

0.49990 

2.00x10^ 

0.847 

0.0513 

Linear 

Elastic 

Specimen 

1.00 

4.30  xKT4 

0.143 

0.144 

0.49950 

0.00100 

0.379 

0.0207 

0.0120 

0.0230 

0.00153 

0.0717 

0.0719 

0.49900 

0.00200 

0.268 

0.0163 

1.44x1  O'4 

0.00717 

0.00736 

0.49000 

0.0201 

0.0858 

0.00520 

1.45x10^ 

0.00358 

0.00378 

0.48000 

0.0405 

0.0615 

0.0121 

0.00373 

These  and  other  relations  between  the  moduli  and  Poisson’s  ratio  may  be  found  in  most 
textbooks  on  elasticity  or  elastic  wave  propagation.  It  is  clear  from  table  1  that  the  aluminum 
pressure  bars  are  substantially  stiffer  than  the  specimen  materials.  This  is  particularly  true  in 
shear.  The  shear  modulus  of  the  specimen  materials  is  4-5  orders  of  magnitude  smaller  than  the 
shear  modulus  of  the  bars. 

The  dimensionless  parameter  G/K,  the  ratio  of  the  shear  to  the  bulk  modulus,  is  a  measure  of  the 
incompressibility  of  the  material.  It  is  a  strictly  decreasing  function  of  Poisson’s  ratio  v  and 
approaches  zero  as  v  and  approaches  Vi  from  below,  as  follows  from  equation  2\.  Materials  for 
which  G/K <  0.0 1  are  typically  referred  to  as  nearly  incompressible .  From  equation  2i,  we  find 
that  this  condition  is  equivalent  to  v>  0.495  4  Note  that  for  the  specimen,  all  cases  but  the  last 
two  in  table  1  satisfy  this  condition.  For  nearly  incompressible  materials,  volumetric  and  shear 
strains  of  the  same  magnitude  will  generate  a  pressure  that  is  several  orders  of  magnitude  larger 
than  the  shear  stress;  conversely,  pressure  and  shear  stresses  of  the  same  magnitude  will  generate 
a  volumetric  strain  that  is  several  orders  of  magnitude  smaller  than  the  shear  strain. 

2.4  Wave  Speeds  and  Travel  Times 

Table  1  also  lists  three  different  linear  elastic  wave  speeds  for  each  material:  the  longitudinal 
wave  speed5  cE,  the  shear  wave  speed  cg,  and  the  “bar  wave  speed”  ce  (see  section  5.1)  The 
subscript  indicates  which  elastic  modulus  is  used  to  compute  the  wave  speed  from  the  relation 
cM  =  -sjM  /  p  ,  where  M  denotes  one  of  the  moduli  L,  G,  or  E.  The  table  does  not  list  the 

Rayleigh  wave  speeds,  which  govern  the  propagation  of  surface  waves.  For  all  the  specimen 
materials  considered  here,  the  Rayleigh  wave  speed  is  about  95%  of  the  corresponding  shear 
wave  speed;  for  the  aluminum  bars,  it  is  about  93%  of  the  shear  wave  speed.6 

The  mechanical  impedances  corresponding  to  the  longitudinal  modulus  and  the  Young’s 
modulus  are  pcL  and  pcE,  respectively.  The  ratio  of  specimen  impedance  to  bar  impedance 
governs  the  transmission  and  reflection  of  the  loading  wave.  Table  1  lists  the  values  of  this  ratio. 
In  the  one-dimensional  theory  of  wave  propagation  the  different  cross-sectional  areas  of  the 
specimen  and  the  bar  would  also  be  taken  into  account,  and  in  such  analyses  the  impedances 
based  on  the  “bar  speed”  (i.e.,  pcE)  are  typically  used.  However,  at  early  times  after  the  arrival 
of  the  loading  wave,  the  specimen  deformation  is  far  from  one-dimensional  and  the  longitudinal 
wave  speed  cE  is  more  relevant,  particularly  in  the  interior  of  the  specimen  and  the  bars. 

Table  2  summarizes  the  geometric  parameters  of  the  problem  and  also  gives  the  corresponding 
wave  travel  times  based  on  the  wave  speeds  in  table  1.  The  travel  times  (in  ps)  are  listed  in  two 
groups:  the  time  for  a  wave  to  travel  the  length  of  the  bar  or  the  specimen,  and  the  time  for  a 


4  Some  authors  use  the  less  restrictive  condition  v  >  0.49  as  the  criterion  for  near  incompressibility;  this  is  equivalent  to  the 
condition  G/K  <  0.02. 

5  This  is  also  referred  to  as  the  dilatational  wave  speed. 

6  See  table  7.5.1  in  Eringen  and  Suhubi  (72). 
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wave  to  travel  the  radius  of  the  bar  or  the  specimen.  Of  course,  these  travel  times  depend  on  the 
type  of  wave  under  consideration.  Travel  times  for  the  lengths  are  given  for  the  longitudinal, 
shear,  and  bar  waves.  Travel  times  for  the  radius  are  given  for  the  longitudinal  and  shear  waves 
only,  since  the  bar  wave  speed  is  not  relevant  in  this  case. 


Table  2.  Geometric  parameters  and  wave  travel  times. 


Material 

Length 

(mm) 

Travel  Times 

(ps) 

Radius 

(mm) 

Travel  Times 
(ps) 

Based 

on  cL 

Based 

on  cE 

Based 
on  cG 

Based 

on  cL 

Based 

oncG 

Incident  Bar 

768.0 

125.7 

153.0 

249.6 

12.8 

2.10 

4.16 

768.0 

125.7 

153.0 

249.6 

9.5 

1.56 

3.09 

Transmission  Bar 

256.0 

41.9 

51.0 

83.2 

12.8 

2.10 

4.16 

256.0 

41.9 

51.0 

83.2 

9.5 

1.56 

3.09 

Ballistic  Gelatin  Specimen 

1.45 

0.72 

93.60 

162.1 

6.35 

3.17 

710.0 

Linear  Elastic  Specimen 
v  =  0.49999 

1.45 

0.54 

69.93 

121.1 

6.35 

2.37 

530.4 

v  =  0.49990 

1.71 

121.1 

7.50 

530.4 

v  =  0.49950 

3.83 

121.1 

16.76 

530.3 

v  =  0.49900 

5.41 

121.1 

23.69 

530.2 

v  =  0.49000 

16.91 

120.7 

74.02 

528.6 

v  =  0.4800 

23.59 

120.3 

103.3 

526.6 

The  wave  speeds  listed  in  tables  1  and  2  for  the  Mooney-Rivlin  model  for  the  specimen  are  only 
approximate  since  they  are  based  on  the  linear  elastic  approximation  to  a  nonlinear  constitutive 
model.  The  wave  speeds,  and  hence  the  impedance  ratios  and  travel  times,  will  change  with  the 
deformation.  Also,  since  the  specimen  undergoes  large  deformations  (for  both  the  nonlinear  and 
linear  elastic  models),  and  since  the  travel  times  are  based  on  the  undeformed  lengths,  these 
times  are  only  approximate  for  the  deformed  specimen.  Nevertheless,  the  tabulated  values 
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provide  an  order  of  magnitude  estimate  for  the  wave  speeds  and  travel  times  in  the  specimen. 
These  estimates  for  the  wave  travel  times  are  useful  for  determining  whether  certain  features  of 
the  specimen  response  (to  be  discussed  later)  correlate  with  the  arrival  of  certain  waves. 

2.5  Measurement  of  Gap  Size 

Many  of  the  figures  in  this  report  give  (for  one  or  both  interfaces)  the  gap  size  as  a  function  of 
time  on  the  centerline  or  the  gap  size  as  a  function  of  radius  at  fixed  times.  In  this  section,  we 
specify  the  convention  used  in  measuring  these  gaps.  We  begin  by  restricting  attention  to  the 
centerline. 

Material  points  initially  on  the  centerline  must  remain  on  the  centerline,  that  is,  r  =  0  if  and  only 
if  R  =  0.  Hence,  the  radial  displacement  of  points  on  the  centerline  is  zero.  The  axial 
displacements  of  the  specimen  and  the  incident  bar  at  their  interface  (Z  =  0)  on  the  centerline  are 
denoted  by  m.s-ib  and  um,  respectively;  similarly,  zs-m  and  z.ui  denote  the  corresponding  deformed 
axial  coordinates.  Since  the  specimen  and  the  bar  are  initially  in  contact,  they  will  continue  to  be 
in  contact  at  some  later  time  t  if  and  only  if  their  axial  displacements  coincide  at  time  t  or, 
equivalently,  their  deformed  axial  coordinates  coincide  at  time  t.  If  the  specimen  and  the 
incident  bar  separate,  that  is,  if  a  gap  forms  between  them,  then  the  size  of  this  gap,  AS-ib,  is 
given  by  the  difference  in  their  displacements  or,  equivalently,  by  the  difference  in  their 
deformed  axial  coordinates: 


As-IB  -  Ws-IB  -  U IB  -  Zs-IB  -  ZlB  • 


(3) 


Similarly,  the  axial  displacements  of  the  transmission  bar  and  the  specimen  at  their  interface 
(Z  =  Ls  =  1.45)  on  the  centerline  are  denoted  by  wT b  and  «s-tb,  respectively;  and  zTb  and  zS-tb 
denote  the  corresponding  deformed  axial  coordinates.  If  the  specimen  and  the  transmission  bar 
separate,  then  the  size  of  the  gap  between  them,  As-tb,  is  given  by  the  difference  in  their 
displacements  or,  equivalently,  by  the  difference  in  their  deformed  axial  coordinates: 

As-TB  =  Wtb  —  Ws-TB  =  ZTB  —  ZS-TB  •  (4) 

As  will  be  seen  in  subsequent  figures,  when  gaps  do  form  it  is  the  result  of  the  specimen  moving 
away  from  the  bars,  so  that  for  the  coordinate  system  in  figure  1,  ws-ib  >  «ib  for  a  gap  at  the 
incident  bar,  and  ws-tb  <  «tb  for  a  gap  at  the  transmission  bar.  The  definitions  of  the  gap  sizes 
given  above  were  chosen  so  that  the  sign  of  the  gap  would  be  positive  in  these  cases;  of  course,  a 
value  of  zero  for  the  gap  implies  that  the  specimen  and  bar  are  in  contact. 

If,  as  above,  we  restrict  attention  to  the  centerline,  then  the  gap  sizes  are  a  function  of  time  only. 
The  majority  of  our  plots  of  gap  size  are  for  this  case.  However,  we  have  also  included  a  few 
plots  of  gap  size  versus  radial  location  at  selected  times.  For  this  case  it  makes  a  difference 
whether  the  axial  displacements  are  referred  to  points  on  the  specimen  and  bar  with  the  same 
initial  radius  R  or  to  points  on  the  specimen  and  bar  with  the  same  deformed  radius  r.  The  latter 
case,  that  is,  the  axial  separation  of  the  specimen  and  bar  along  a  line  of  constant  radius  in  the 
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deformed  state,  is  of  more  interest  and  is  used  here.  With  this  caveat,  the  definitions  above  still 
apply.7 

2.6  Stress,  Stretch  and  Strain 

When  the  specimen  remains  in  contact  with  the  bars,  the  deformed  thickness  (or  length)  of  the 
specimen,  £s,  is  determined  either  from  the  current  positions  of  its  faces  or  from  its  initial  length 
Ls  and  the  current  displacements  of  its  faces: 

£s  =  Ztb  -  Sib  =  Ls  +  Wtb  -  z<ib  .  (5) 

These  variables  may  be  measured  at  any  radial  location  since  the  faces  of  the  bars  remain  very 
nearly  planar  when  the  specimen  is  relatively  soft,  as  in  the  cases  considered  here.  At  radial 
locations  where  there  are  substantial  gaps  between  the  specimen  and  the  bars,  equation  5  will 
overestimate  the  deformed  thickness. 


The  term  stretch  denotes  the  local  ratio  of  deformed  to  undeformed  length,  which,  of  course, 
varies  with  direction.  The  stretch  in  all  directions  is  unity  in  the  undeformed  state.  A  stretch 
greater  than  1  corresponds  to  extension  in  the  given  direction;  a  stretch  between  0  and  1 
corresponds  to  compression.  The  stretches  in  the  radial,  hoop  and  axial  directions  are  denoted  by 
Xr,  Xq  and  Xz,  respectively.  They  are  the  diagonal  components  of  the  deformation  gradient  tensor 
F  and  are  given  by  the  relations 


4. =4. 


R 


dZ 


(6) 


We  are  primarily  interested  in  the  stretches  and  strains  in  the  specimen  and,  in  particular,  in  the 
axial  stretch  and  corresponding  axial  strain.  The  nominal  (or  engineering)  strain  e7  and  the 
logarithmic  (or  true)  strain  sz  in  the  axial  direction  are  given  by 

ez  —  1—  Xz ,  Sj  —  —In  Xz .  (7) 

Both  of  these  strain  measures  are  used  in  the  experimental  literature  for  reporting  SHBP  test 
data.  The  radial  and  hoop  strains  are  defined  similarly.  Note  that  these  strain  measures  have 
been  taken  positive  if  the  particular  coordinate  direction  is  in  compression.  In  the  general 
continuum  mechanics  literature,  the  sign  convention  for  strain  is  “positive  in  extension”  and 
“negative  in  compression”.  However,  in  the  SHPB  literature  the  opposite  convention  is  typically 
used  since  the  specimen  is  under  axial  compression;  we  have  followed  that  convention  here. 

The  average  or  mean  value  of  the  axial  stretch  in  the  specimen  is  the  ratio  of  the  deformed  to 
undeformed  thickness: 


7  The  deformed  axial  coordinates  are  computed  at  element  nodes.  The  specimen  and  bar  nodes  off  the  centerline  do  not 
necessarily  line  up.  For  a  specimen  node  with  deformed  radius  r,  we  used  linear  interpolation  between  the  two  incident  bar  nodes 
whose  deformed  radius  bounded  r  in  order  to  determine  zIB(r);  similarly  for  the  transmission  bar. 
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(8) 


_  WTB  MIB 


u 


where  the  relation  on  the  right  follows  from  equation  5.  The  mean  value  of  the  nominal  axial 
strain8  is  given  by 


e  =1-A  =  1-^-  =  Um  “™ 


(9) 


Note  that  we  use  the  same  notation  for  local  and  mean  values;  it  should  be  clear  from  the  context 
which  meaning  is  intended. 

We  use  a  superposed  dot  to  denote  the  material  time  derivative.  It  follows  from  equation  9  that 
the  mean  value  of  the  nominal  axial  strain  rate  is  given  by 


(10) 


where  v®  and  vtb  denote  the  axial  component  of  velocity  of  the  incident  and  transmission  bars  at 
the  specimen  interfaces.  The  approximation  on  the  right  follows  from  the  fact  that  for  soft 
specimens,  vib  will  generally  be  much  larger  than  vtb-9  This  approximate  relation  yields  an  even 
cruder  approximation  for  the  “strain  acceleration”:  e_  =  -A.  «  vm  /  L§ .  For  the  imposed  velocity 

history  at  the  far  end  of  the  incident  bar,  as  described  in  equation  1,  v®  eventually  oscillates 
around  a  mean  value  of  roughly  twice  v0  (see  section  5.3). 

Note  that  it  is  the  mean  values  of  axial  strain  and  axial  strain  rate  that  are  actually  measured10  in 
an  SHPB  test.  Similarly,  the  values  of  the  axial  stretch  given  in  our  figures  are  the  mean  values 
computed  from  equation  8;  consequently,  the  corresponding  nominal  axial  strains  cited  in 
subsequent  discussions  are  the  mean  values  determined  from  equation  9.  These  mean  values 
may  differ  substantially  from  the  local  values  if  the  specimen  deformation  is  highly  non-uniform. 
This  will  certainly  be  the  case  during  the  “ring-up”  period  when  wave  propagation  dominates. 
Also,  at  radial  locations  where  there  are  substantial  gaps  between  the  specimen  and  the  bars,  the 
relation  8  will  overestimate  the  axial  stretch,  and  consequently,  the  relation  9  will  underestimate 
the  axial  strain.  Likewise,  the  relations  for  the  mean  value  of  the  nominal  axial  strain  rate  in 
equation  10  will  not  be  valid  when  a  gap  forms. 

The  Cauchy  stress  tensor  is  denoted  by  a.  Its  components,  which  give  force  per  unit  deformed 
area,  are  often  referred  to  as  true  stresses.  We  follow  the  sign  convention  in  the  general 


8  A  mean  value  of  the  logarithmic  axial  strain  is  provided  by  equation  72  with  Xz  given  by  equation  8. 

0  The  displacement  histories  at  the  S-IB  and  S-TB  interfaces  in  figures  7  and  8  are  consistent  with  this  assertion. 

10  More  precisely,  they  are  inferred  from  strain  gage  measurements  in  the  pressure  bars. 
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continuum  mechanics  literature,  that  is,  the  stress  components  are  taken  positive  in  tension  and 
negative  in  compression,11  whereas  the  pressure  p  is  taken  to  be  positive  in  compression: 


1  1/  X 

Z1  =  -  3 tr  =  -  3  Kr  +  °e0  +  • 


(11) 


Here,  “tr”  denotes  the  trace,  and  a)T,  aee,  and  ozz  are  the  radial,  hoop  and  axial  components 
of  a.  Since  axisymmetric  deformations  are  imposed  and  the  constitutive  models  for  all  materials 
are  isotropic, 

CT0J —  OjQ —  CJ;-0 —  CT@; —  0.  (12) 

However,  the  shear  stress  component  orz  =  ozr  need  not  be  zero  in  general;  in  particular,  orz  will 
be  nonzero  in  the  specimen  whenever  there  is  bulging. 

Since  only  isotropic  elastic  constitutive  relations  are  used  here,  the  principal  axes  of  stress  and 
strain  coincide.  The  9  (i.e.,  hoop)  direction  is  always  a  principal  axis  of  stress  and  strain  for 
axisymmetric  deformations,  but  the  r  and  z  (i.e.,  radial  and  axial)  directions  will  be  principal 
axes  if  and  only  if  the  shear  stress  orz  is  zero.  In  this  case,  the  radial,  hoop  and  axial  stretches  are 
the  principal  stretches,  and  the  Jacobian12  of  the  deformation,  J,  is  given  by 

J  =  A,-AqAz.  (13) 

There  are  two  situations  where  the  condition  orz  =  0,  the  relation  13,  and  the  relations 

orr  =  cj00 ,  er  =  sq  ,  Ar  —  Aq  —  r- —  ,  (14) 


can  be  guaranteed  to  hold  exactly.  One  of  these  occurs  when  the  deformation  is  uniform,  which 
is  the  desired  state  in  an  SHPB  test.  The  other  occurs  at  any  point  on  the  centerline  of  a  solid 
specimen,  regardless  of  whether  or  not  the  deformation  is  uniform.13  Furthermore,  for  the  nearly 
incompressible  specimen  materials  considered  here,  J  will  generally  be  close  to  one  even  for 
large  deformations,  in  which  case  equation  143  yields  the  approximation 


Ar  - 


•K 


(15) 


11  This  is  also  the  sign  convention  used  in  LS-DYNA.  However,  in  the  SHPB  literature  it  is  more  common  to  take  stress 
positive  in  compression,  since  the  specimen  is  under  axial  compression  during  the  useful  part  of  the  test. 

-j  O 

This  is  the  determinant  of  the  deformation  gradient  F  and  represents  the  local  ratio  of  deformed  to  undeformed  volume. 
13  That  arr  =  aee  and  ara=  0  on  the  centerline  follow  from  the  radial  and  axial  components,  respectively,  of  the  momentum 
balance  equation,  after  multiplying  by  r  and  taking  the  limit  as  r  approaches  zero.  That  Ar  =  Ag  on  the  centerline  follows  from 

equations  6,  since  r/R  approaches  drl dR  as  R  and  r  approach  zero.  Also,  both  drj dZ  and  dz/ SR  must  be  zero  on  the 
centerline  to  prohibit  singularities. 
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on  the  centerline.  When  the  deformation  in  the  specimen  is  non-uniform,  the  relations  13-15  will 
hold  approximately  at  points  near  the  centerline,  and  they  may  even  hold  approximately  at  points 
substantially  off  the  centerline. 


3.  The  Mooney-Rivlin  Constitutive  Model  for  the  Specimen 


3.1  Incompressible  and  Compressible  Versions 

The  classical  Mooney-Rivlin  model  is  a  nonlinear  elastic  constitutive  model  for  isotropic, 
incompressible  solids.  It  is  one  of  the  earliest  models  of  this  type,  and  was  originally  applied  to 
the  response  of  rubber  under  large  deformations;  see  the  original  papers  by  Mooney  (13)  and 
Rivlin  (14)  and  the  book  by  Treloar  (15).  Because  of  its  simplicity,  the  model  is  widely  used  in 
theoretical  studies  and  is  often  applied  to  nearly  incompressible  materials  other  than  rubber.  The 
Mooney-Rivlin  model  is  discussed  in  many  textbooks  on  nonlinear  elasticity  and  continuum 
mechanics;  see  (16-19).  A  description  of  the  model  is  provided  in  appendix  A-2. 

The  incompressibility  assumption  is,  of  course,  an  idealization.  For  quasi-static  problems  with  at 
least  one  free  surface,  the  assumption  of  incompressibility  typically  yields  reasonable  results  if 
the  bulk  modulus  is  at  least  two  orders  of  magnitude  larger  than  the  shear  modulus.  However, 
the  incompressibility  constraint  is  inconsistent  with  the  propagation  of  longitudinal  waves,  so  for 
dynamic  problems  involving  nearly  incompressible  materials,  some  degree  of  compressibility 
must  be  added  to  the  constitutive  model. 

The  simulations  described  in  sections  6  and  7  and  in  appendices  B  and  C  used  LS-DYNA’s 
compressible  version  of  the  Mooney-Rivlin  model  for  the  specimen;  this  is  Material  Type  27  in 
the  manuals  (8,  9).  The  model  has  three  material  constants  (aside  from  the  density):  the  elastic 
moduli14  Ai  and  A 2,  and  the  Poisson’s  ratio  v.  The  two  moduli  characterize  the  response  to 
volume-preserving  deformations  (just  as  in  the  incompressible  version),  whereas  the  nonlinear 
pressure-volume  relation  involves  all  three  constants.  We  refer  to  this  as  the  “compressible 
Mooney-Rivlin  model”;  it  is  discussed  in  more  detail  in  appendix  A-3.  Verification  of  the  model 
implementation  is  discussed  in  appendix  A-4.  For  the  purposes  of  the  present  discussion,  we 
simply  note  that  the  shear  modulus  G  in  the  small  strain,  linear  elastic  approximation  to  the 
model  is  related  to  the  moduli  A\  and  A2  by 

G  =  2(Al+A1).  (16) 

Since  the  linear  elastic  bulk  modulus  K  is  given  in  terms  of  the  shear  modulus  and  Poisson’s 
ratio  by  equation  2],  v  can  be  adjusted  so  as  to  yield  the  bulk  modulus  for  the  material  in 
question,  once  Ai  and  A2  (and  hence  G)  have  been  determined.  For  the  nearly  incompressible 
materials  of  interest  here,  this  means  choosing  v  very  close  to  V2. 

14  The  LS-DYNA  manuals  (8,  9)  use  the  symbols  A  and  B  for  the  constants  A{  and  A2,  respectively. 
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3.2  Calibration  of  the  Mooney-Rivlin  Model  for  Ballistic  Gelatin 

Ballistic  gelatin  is  a  prime  example  of  a  soft,  nearly  incompressible  material.  Gelatin  powder  is 
an  aggregate  of  large  protein  molecules  (collagen)  of  various  sizes,  extracted  from  animal  tissues 
(20).  In  this  report,  the  term  gelatin  refers  to  the  solution  of  this  powder  in  water.  The  common 
formulations  of  ballistic  gelatin  are  20%  and  10%  powder  by  mass  (80%  and  90%  water  by 
mass).  Both  formulations  are  used  as  tissue  surrogates  in  impact  and  penetration  tests. 

LS-DYNA’s  compressible  Mooney-Rivlin  model  was  calibrated  to  yield  a  response  typical  of 
20%  ballistic  gelatin.  The  parameters  A\  and  Ai_  listed  in  table  3  were  chosen  to  give  rough 
agreement  with  quasi-static  uniaxial  compression  data  on  ballistic  gelatin  obtained  in  an 
experimental  study  by  Moy  et  al.  (5).  More  precisely,  these  values  give  a  stress-strain  curve  that 
lies  slightly  above  their  curve  for  a  strain  rate  of  1/s.15  This  represents  a  rough  extrapolation  of 
their  quasi-static  data  (strain  rates  of  0.001/s,  0.01/s,  1/s)  to  strain  rates  on  the  order  of  2500/s.16 
By  equation  16,  the  corresponding  small  strain  shear  modulus  is  G  =  80  kPa.  Substantially  lower 
values  for  these  moduli  would  be  expected  for  10%  gelatin. 


Table  3.  Parameters  for  the  compressible  Mooney-Rivlin  Model. 


p  (g/cm3) 

A,  (kPa) 

A2  (kPa) 

v 

1.00 

12.0 

28.0 

0.49999 

Since  gelatin  is  mostly  water,  one  would  expect  its  density  and  (initial)  bulk  modulus  to  be  close 
to  that  of  water,  namely,  1  g/cm  and  2.3  GPa,  respectively.  The  room  temperature  densities  for 
20%  gelatin  reported  in  Aihaiti  and  Hemley  (10)  and  Winter  and  Shifler  (20)  are  about 
1.00  g/cm3  and  1.06  g/cm3,  respectively.17  We  used  the  value  1.00  g/cm3.  However,  the 
estimated  bulk  modulus  for  20%  ballistic  gelatin  in  Aihaiti  and  Hemley  (10)  is  more  than  twice 
that  of  water.18  We  chose  a  Poisson’s  ratio  of  v  =  0.49999,  which  together  with  the  shear 
modulus  of  80  kPa  (calculated  above)  yields  a  bulk  modulus  of  K  =  4.0  GPa.  Note  that  this  is 
more  than  five  orders  of  magnitude  larger  than  the  shear  modulus. 


15  The  corresponding  stress-stretch  curve  is  plotted  in  figure  A-l  in  our  appendix  A-4.  As  discussed  their,  these  curves  are 
largely  insensitive  the  value  of  v,  provided  that  it  is  close  to  Vi. 

16  This  extrapolated  curve  lies  well  below  their  stress-strain  curve  at  a  strain  rate  of  2500/s,  obtained  in  an  SHPB  test  on  an 
annular  specimen.  We  believe  that  there  are  large  inertial  effects  in  that  test  and,  in  particular,  that  the  stress  state  is  not  uniaxial 
and  so  not  comparable  with  the  quasi-static  data.  The  quasi-static  data  referred  to  here  was  provided  by  Tusit  Weerasooriya  of 
ARL.  This  data  is  not  included  in  (5),  except  for  the  upper  stress-strain  curve  in  figure  6(a)  in  that  report,  which  is  described 
there  as  corresponding  to  a  strain  rate  of  1/s  but  which  appears  to  be  the  0.001/s  curve  in  the  data  provided  to  us. 

17  The  value  1 .00  g/cm3  is  based  on  a  slight  extrapolation  of  the  data  given  in  Aihaiti  and  Hemley  (10):  the  authors  attribute 
this  data  to  Dana  Dattlebaum  at  Los  Alamos  National  Laboratory  (LANL). 

18  The  value  reported  in  Aihaiti  and  Hemley  (10)  varied  from  about  4.6  to  5.2  GPa  at  room  temperature,  depending  on  the 
functional  form  of  the  equations  used  to  fit  the  data.  However,  it  appears  that  even  the  lowest  of  these  values  may  be  too  high. 
After  are  computational  study  was  completed,  we  were  informed  that  the  data  in  Aihaiti  and  Hemley  (10)  is  inconsistent  with  the 
shock  Hugoniot  data  from  LANL,  and  that  the  latter  implies  an  initial  bulk  modulus  much  closer  to  that  of  water  (21). 
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We  emphasize  that  our  intention  in  this  study  was  neither  to  develop  a  constitutive  model  for 
ballistic  gelatin  nor  to  simulate  a  particular  SHPB  experiment  on  gelatin.  Ballistic  gelatin  is 
known  to  exhibit  viscoelastic  behavior.  The  rate-dependence  of  the  large  strain  compression 
data  in  (5)  and  elsewhere  is  consistent  with  this,  as  is  the  small  strain  rheometer  data  in  Juliano  et 
al.  (22).  An  accurate  model  for  ballistic  gelatin  would  require  the  addition  of  a  viscoelastic 
component  to  the  constitutive  model.  Consequently,  numerical  simulations  capable  of 
duplicating  high-rate  SHPB  experiments  would  likely  require  a  properly  calibrated  nonlinear 
viscoelastic  model.  The  use  of  a  purely  elastic  model  was  motivated  by  our  original  objective  of 
quantifying  inertial  effects  in  SHPB  tests.  Since  there  are  no  strain  rate  effects  in  an  elastic 
model,  any  differences  observed  between  quasi-static  and  high-rate  uniaxial  compression  tests 
must  be  due  to  inertial  effects  in  the  latter. 


4.  Computational  Parameters 


4.1  Geometry  and  Meshing 

The  meshes  for  the  specimen  and  pressure  bars  consisted  of  4-node  quadrilateral  elements  used 
in  2-D  axisymmetric  mode.  Single-point  integration  was  used  throughout.  Three  different 
meshes  were  used  for  the  specimen,  as  shown  in  figure  2.  Since  the  axial  strain  in  the  specimen 
is  approximately  twice  the  radial  strain,  it  was  beneficial  to  use  specimen  elements  whose  initial 
length  in  the  axial  direction  was  twice  that  in  the  radial  direction,  so  that  the  element  became 
more  nearly  square  under  severe  axial  compression.  The  default  or  baseline  mesh  is  shown  in 
figure  2a  and  referred  to  as  the  “25x50  pm  mesh”:  all  specimen  elements  had  a  radial  edge 
length  of  25  pm  and  an  axial  edge  length  of  50  pm.  This  specimen  mesh  was  used  in  all 
simulations  with  the  exception  of  the  mesh  sensitivity  study  discussed  in  appendix  B-l,  where 
the  results  for  the  baseline  mesh  are  compared  with  the  results  for  a  coarser  mesh  (50x100  pm  in 
figure  2b)  and  a  finer  mesh  (12.3x12.5  pm  in  figure  2c) 
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Figure  2.  Meshes  for  the  Incident  Bar  (IB),  Transmission  Bar  (TB),  and  Specimen  (S  ).  The 
baseline  mesh  for  the  specimen  is  shown  in  case  (a). 

All  incident  bar  and  transmission  bar  elements  were  squares  with  a  1 00- pm  edge  length. 
Differences  between  bar  and  specimen  element  sizes  in  the  axial  direction  could  result  in  a 
computational  impedance  mismatch  if  the  specimen  and  bars  had  similar  mechanical 
impedances.  This  was  not  an  issue  here  since  the  specimen  and  bar  impedances  were  already 
quite  dissimilar  (see  table  1).  Also,  since  the  faces  of  the  stiffer  bars  are  expected  to  remain 
nearly  flat,  elements  with  larger  radial  edge  length  could  be  tolerated.  The  element  size  in  the 
bars  does  need  to  be  small  enough  to  capture  the  essential  features  of  the  wave  propagation  (see 
section  5).  The  100- pm  edge  length  seemed  to  be  sufficient  for  this. 

4.2  Boundary  Conditions 

The  interaction  of  the  specimen  and  the  bars  at  the  S-IB  interface  and  the  S-TB  interface  was 
governed  by  the  LS-DYNA  segment-based  penalty  algorithm  CONTACT_2D  with  the  option 
AUTOMATIC_SURFACE_TO_SURFACE  (8,  9).  This  algorithm  requires  three  input 
parameters:  SFACT,  a  scale  factor  for  the  penalty  force  stiffness;  VDC,  a  viscous  damping 
coefficient  (in  percent  of  the  level  for  critical  damping);  and  a  Coulomb  friction  coefficient, 
which  was  set  to  zero  for  frictionless  contact  (see  section  2.2).  The  LS-DYNA  default  values 
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SFAC  =  1  and  VDC  =  10  were  used  for  all  simulations  with  the  exception  of  the  contact 
algorithm  parameter  sensitivity  study  in  appendix  B-2. 

As  discussed  in  section  2.2,  the  boundary  condition  at  the  impact  end  of  the  incident  bar  is  a 
prescribed  axial  velocity  history.  No  radial  constraint  was  placed  on  the  nodes  there. 

At  the  far  end  of  the  transmission  bar,  we  used  the  LS-DYNA  non-reflecting  boundary  condition 
BOUND ARY_NON_REFLECTING_2D  (5,  9)  to  simulate  a  semi-infinite  transmission  bar. 

This  algorithm  imposes  normal  and  shearing  stresses  given  by 

^normal  P  GY  normal  ’  ^shear  P  U?  ^tangential  '  , 

Here  cl  and  cq  are  the  longitudinal  and  shear  wave  speeds,  respectively,  of  the  transmitting 
material  (in  our  case  aluminum),  and  vnormai  and  vtangentiai  are  the  particle  velocities  in  the  normal 
and  tangential  direction,  respectively.  We  conducted  some  preliminary  tests  on  this  boundary 
condition  and  found  that  it  performed  extremely  well  in  one-dimensional  uniaxial  strain,  and 
reasonably  well  for  waves  of  the  type  generated  in  our  simulations.  Once  the  loading  wave 
travels  through  the  specimen  and  reaches  the  transmission  bar,  it  takes  about  42  ps  for  a 
longitudinal  wave  and  51  ps  for  the  bar  wave  to  travel  the  length  of  the  transmission  bar  (see 
table  2),  and  then  another  42-51  ps  for  any  slight  reflections  from  the  supposedly  non-reflecting 
boundary  to  arrive  back  at  the  specimen.  Thus  any  influence  that  slight  reflections  from  the  far 
end  of  the  transmission  bar  might  have  on  the  interactions  at  the  S-TB  interface  could  not 
possibly  be  felt  earlier  than  84-102  ps  after  the  loading  wave  first  reached  that  interface. 

4.3  Miscellaneous 

The  default  values,  1.5  and  0.06,  were  used  for  the  quadratic  and  linear  artificial  viscosity 
coefficients,  respectively.  The  time  step  factor  of  safety,  TSSFAC,  was  set  to  0.1,  meaning  that 
the  time  step  computed  using  the  Courant  condition  is  then  multiplied  by  0.1. 

LS-DYNA’s  Orthotropic  Elastic  model  (Material  Model  2  in  LS-DYNA  [§]  and  Hallquist  [9]) 
was  used  for  the  pressure  bars.  The  material  constants  for  this  orthotropic,  linear  elastic  model 
were  chosen  in  such  a  way  as  to  yield  an  isotropic,  linear  elastic  model.  Values  of  the  material 
constants  appropriate  for  aluminum  were  used  in  all  the  simulations  and  are  given  in  table  1. 


5.  Wave  Propagation  in  the  Pressure  Bars 


5.1  General  Remarks 

The  loading  wave  in  the  incident  bar  generates  an  axial  stress  on  the  order  of  pceVo,  where  Vo  is 
the  plateau  of  the  imposed  velocity  at  the  impact  end,  as  described  in  equation  1.  For  the  value 
vq  =  1.8125  m/s  used  in  our  simulations  this  gives  an  axial  stress  on  the  order  of  25  MPa  for 
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aluminum,  which  is  well  below  the  dynamic  yield  stress.19  This  justifies  the  use  of  a  linear 
elastic  model  for  the  pressure  bars  in  our  simulations.  Linear  elastic  deformation  of  the  bars  is 
also  a  requirement  in  actual  SHPB  tests,  since  the  analysis  of  test  data  is  based  on  linear  elastic 
wave  propagation  theory  ( 1 ,  2). 

A  good  reference  for  linear  elastic  wave  propagation  in  bars  is  Chapters  2  and  8  in  Graff  (24); 
see  also  Eringen  and  Suhubi  (12).  The  one-dimensional  theory  of  wave  propagation  in  long  bars 
predicts  that  a  longitudinal  elastic  wave  travels  undistorted  at  the  “bar  speed”  ce,  that  is,  at  the 
speed  governed  by  the  Young’s  modulus  E.  The  more  accurate  results  of  the  three-dimensional 
theory  are  much  more  complex,  particularly  near  the  loading  end  where  the  wave  becomes 
highly  distorted  with  large  oscillations.  These  oscillations  decay  somewhat  as  the  wave 
progresses  down  the  bar,  but  they  do  not  go  away  entirely.  Near  the  centerline  of  the  bar  the 
wave  initially  propagates  at  the  longitudinal  elastic  wave  speed  cl .  But  release  waves  from  the 
stress-free  lateral  surface  interact  with  the  longitudinal  wave,  so  that  by  the  time  the  wave  has 
traveled  10-20-bar  diameters  only  a  trace  of  the  initial  longitudinal  wave  remains.  Most  of  the 
energy  is  contained  in  a  wave  moving  at  the  slower  bar  speed  ce ■  Due  to  wave  dispersion,  the 
rise  time  of  the  loading  pulse  increases  as  it  moves  down  the  bar  and  the  axial  stress,  strain  and 
particle  velocity  oscillate  about  the  value  predicted  by  the  one-dimensional  theory. 

According  to  table  2,  for  the  material  properties  and  incident  bar  length  considered  here,  the 
longitudinal  and  bar  waves  should  take  about  126  and  153  ps,  respectively,  to  arrive  at  the 
specimen.  The  average  of  these  arrival  times  is  about  140  ps.  As  will  be  seen  in  section  5.3,  the 
velocity  histories  in  the  incident  bar  in  the  simulations  are  consistent  with  the  preceding 
observations. 

5.2  Considerations  Governing  the  Bar  Lengths 

In  an  SHPB  experiment  the  length  of  the  striker  bar  determines  the  duration  of  the  loading  pulse, 
which  in  turn  determines  the  duration  of  the  test.  The  length  of  the  incident  bar  is  determined  by 
the  requirement  that  the  incident  and  reflected  pulses  measured  by  a  strain  gage  mounted 
midway  down  the  bar  must  not  overlap.  As  a  result,  incident  bar  lengths  of  several  meters  are 
commonly  used. 

In  order  to  reduce  the  size  of  the  output  files  and  particularly  the  run  times  of  the  simulations,  we 
did  not  wish  to  model  the  full  length  of  a  real  incident  bar.  Indeed,  there  was  no  need  to  do  so, 
as  the  prescribed  velocity  history  (equation  1)  is  not  a  pulse  but  instead  has  essentially  infinite 
duration,  and  there  is  no  strain  gage  to  consider.  Nevertheless,  two  considerations  govern  the 
minimum  length  required  for  the  incident  bar  in  our  numerical  simulations.  First,  a  bar  of  length 
at  least  20  diameters  is  desirable  in  order  to  simulate  the  effects  of  dispersion  on  the  loading 
wave.  Second,  in  the  simulations  the  total  strain  that  can  be  imposed  on  the  specimen  is 
proportional  to  the  length  of  incident  bar.  To  see  this,  note  that  the  compressive  incident  wave 
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The  value  of  the  Hugoniot  elastic  limit  for  7075-T6  aluminum  listed  in  Steinberg  (23)  is  420  MPa. 
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reflects  from  the  lower  impedance  specimen  as  a  tensile  wave,  travels  back  down  the  bar,  reflects 
from  the  loading  end  as  tensile  wave,20  and  eventually  arrives  back  at  the  specimen,  after  which 
time  the  simulation  is  no  longer  useful.  The  time  for  this  round  trip  is  thus  twice  the  time  to 
travel  one  length  of  the  bar,  and  this  time  window  limits  the  total  strain  that  can  be  imposed  on 
the  specimen. 

We  chose  an  incident  bar  length  of  Lib  =  768  mm,  which  gives  a  length-to-diameter  ratio  of  30.0 
for  the  larger  diameter  bar  and  a  time  window  of  roughly  2x  140  ps  =  280  ps,  which  suffices  for 
large  strains  in  the  specimen.  Since  the  incident  wave  arrives  at  the  specimen  at  about 
t  =  140  ps,  specimen  loading  would  occur  from  about  t  =  140  ps  to  t  =  420  ps.  However,  when 
gaps  formed  (as  they  did  in  most  of  our  simulations),  they  developed  shortly  after  the  arrival  of 
the  incident  wave.  And  since  this  became  the  focus  of  our  study,  we  terminated  most  of  the 
simulations  at  340  ps.  Thus  most  of  our  time  history  plots  are  restricted  to  the  200-ps  time 
window  from  t  =  140  ps  to  t  =  340  ps. 

As  discussed  at  the  end  of  section  4.2,  slight  reflections  from  the  supposedly  non-reflecting 
boundary  at  the  far  end  of  the  transmission  bar  could  not  possibly  be  felt  at  the  S-TB  interface 
earlier  than  84-102  ps  after  the  loading  wave  first  reaches  that  interface.  Using  t  =  140  ps  for 
the  latter  arrival  time,21  we  see  that  the  effects  of  any  deficiencies  in  the  non-reflecting  boundary 
condition  cannot  be  felt  earlier  than  about  t  =  230  ps.  As  will  be  seen  below,  gaps  formed  at  the 
S-TB  interface  well  before  this  time,  although  we  cannot  rule  out  the  possibility  of  some 
influence  at  later  times. 

5.3  Velocity  Histories  in  the  Incident  Bar 

Figures  3-6  are  plots  of  the  time  history  of  the  axial  component  of  the  particle  velocity,  vz,  at 
four  locations  on  the  centerline  of  the  incident  bar  (r  =  R  =  0).  The  velocity  is  scaled  by  the 
plateau  value  (vo  =  1.8125  m/s)  of  the  input  velocity.  All  plots  are  for  the  bar  with  the  larger 
radius  (7?B  =12.8  mm),  although  the  results  for  the  smaller  radius  should  be  similar. 

Figures  3  and  4  give  velocity  histories  at  the  loading  end  (Z  =  -Lib),  1/3  of  the  way  from  the 
loading  end  (Z  =  -0.666  Lib),  and  2/3  of  the  way  from  the  loading  end  (Z  =  -0.333  Lib).  On 
these  two  plots  the  time  is  scaled  by  the  transit  time  of  a  longitudinal  wave  along  the  length  of 
the  bar,  which  is  Lib  /c^1  =  125.7  ps,  where  c^1  is  the  longitudinal  wave  speed  of  aluminum. 

Figure  3  is  for  an  initial  rise  time  tR  of  1  ps,  which  is  too  short  for  the  linear  ramp  to  be 
distinguishable  on  the  time  scale  of  this  plot.  Figure  4  is  for  the  25-ps  initial  rise  time;  in  this 
case  the  linear  ramp  is  clearly  distinguishable.  For  both  cases,  the  velocity  at  the  two  interior 
locations  oscillates  about  the  plateau  value  vq. 


20  The  fixed  velocity  boundary  condition  at  the  loading  end  produces  the  same  effect  on  the  stress  wave  as  a  rigid  boundary 
would,  that  is,  a  tensile  wave  is  reflected  as  a  tensile  wave.  Indeed,  relative  to  an  observer  moving  with  the  imposed  velocity  v0, 
that  boundary  would  appear  rigid. 

^  Note  that  the  travel  time  through  the  specimen  is  less  than  1  ps  for  the  Mooney-Rivlin  model  as  well  as  for  the  linear 
elastic  model  with  v  =  0.49999  (see  table  2). 
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Figure  3.  The  scaled  axial  velocity  vs.  non-dimensional  time  at  three  locations  on  the  centerline  of  the 

incident  bar:  the  loading  end  with  a  l-|as  rise  time,  1/3  of  the  way  from  the  loading  end.  and  2/3 
of  the  way  from  the  loading  end  (blue  curve). 
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Figure  4.  The  scaled  axial  velocity  vs.  non-dimensional  time  at  three  locations  on  the  centerline  of  the 

incident  bar:  the  loading  end  with  a  25-|is  rise  time,  1/3  of  the  way  from  the  loading  end.  and  2/3 
of  the  way  from  the  loading  end  (blue  curve). 
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Figure  5.  The  scaled  axial  velocity  history  on  the  centerline  of  the  incident  bar  at  the  S-IB  interface  for  a 
1-ps  (initial)  rise  time. 

Figures  5  and  6  give  the  velocity  history  on  the  centerline  of  the  incident  bar  at  the  S-IB  interface 
for  the  1-  and  25-|us  initial  rise  times,  respectively.22  In  these  two  figures  the  time  is  given  in  ps 
for  later  comparison  with  other  plots.  Note  that  the  particle  velocity  ramps  up  to  nearly  double 
its  initial  value  plateau  v0,  as  would  be  expected  from  reflection  of  the  wave  at  a  free  surface. 
This  is  consistent  with  the  small  specimen/bar  impedance  ratio  (see  table  1).  This  doubling  of 
the  velocity  on  reflection  is  also  evident  at  later  times  for  the  two  interior  locations  in  figures 
3  and  4.  The  Mooney-Rivlin  model  was  used  for  the  specimen  in  these  simulations,  but  in  view 
of  the  low  impedance  ratios  for  all  specimen  models  considered  here,  only  small  differences  in 
the  particle  velocity  histories  in  the  incident  bar  after  reflection  would  be  expected  for  the  other 
specimen  models. 


22  The  incident  bar  velocity  v.  in  these  figures  was  denoted  by  vIB  in  equation  10. 
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Figure  6.  The  scaled  axial  velocity  history  on  the  centerline  of  the  incident  bar  at  the  S-IB  interface  for  a 
25-|as  (initial)  rise  time. 

A  comparison  of  figures  3  and  5  (1-ps  initial  rise  time)  with  figures  4  and  6  (25 -ps  initial  rise 
time)  shows  that  the  computational  “pulse  shaping”  has  substantially  smoothed  the  loading  wave 
in  the  latter  case.  Also  of  interest  is  the  increase  in  the  final  rise  times  of  these  waves  due  to 
dispersion  in  the  bar.  This  final  rise  time  could  be  defined  as  the  difference  between  the  time  at 
the  first  peak  in  velocity  at  the  S-IB  interface  and  the  arrival  time  of  the  wave  at  the  S-IB 
interface.  However,  precise  arrival  times  are  difficult  to  determine  from  the  plots  in  figures 
5  and  6.  A  close  inspection  of  figure  5,  for  example,  reveals  barely  discernable  peaks  in  the 
oscillatory  toe  of  the  velocity  history  starting  at  about  135  ps,  although  there  is  no  substantial 
increase  in  velocity  until  about  140  ps.  Thus,  for  the  1-ps  initial  rise  time,  most  of  the  final  rise 
in  particle  velocity  at  the  S-IB  interface  occurs  from  140  to  164  ps,  that  is,  over  a  24-ps  time 
interval.  For  the  25-ps  initial  rise  time,  most  of  the  final  rise  in  particle  velocity  at  the  S-IB 
interface  occurs  from  143  ps  to  184  ps,  that  is,  over  a  41-ps  time  interval. 
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Of  the  three  curves  in  figures  3  and  4,  the  blue  curves  correspond  to  a  point  in  the  incident  bar 
closest  to  the  specimen,  namely,  1/3  of  the  bar  length  back  from  the  S-IB  interface.  In  the 
following  discussion,  we  consider  the  velocity  history  at  this  interior  location  at  early  times,  that 
is,  prior  to  the  arrival  of  the  reflected  wave  at  a  non-dimensional  time  of  about  1.6.  This 
represents  the  loading  wave  shortly  before  it  arrives  at  the  S-IB  interface.  For  a  given  initial  rise 
time,  we  compare  this  loading  wave  with  the  velocity  history  at  the  S-IB  interface  over  a  time 
interval  that  includes  the  first  two  major  peaks  in  velocity  at  these  two  locations.  Clearly,  these 
velocity  histories  are  qualitatively  similar.23  In  particular,  for  the  25 -ps  rise  time,  a  close 
inspection  of  the  plots  reveals  that  both  curves  have  three  inflection  points  prior  to  the  first  peak 
in  velocity.  For  the  1-ps  rise  time,  the  situation  is  more  complicated  due  to  the  small  oscillations 
superimposed  on  the  main  curves;  but  if  these  are  neglected,  both  curves  would  appear  to  have 
only  one  inflection  point  between  the  toe  and  the  first  peak.  These  inflection  points  in  the 
velocity  history  correspond  to  extreme  values  in  the  acceleration  of  the  bar,  which  turn  out  to  be 
important  for  interpreting  the  stress  histories  in  the  specimen.  Furthermore,  the  fact  that  these 
features  are  also  exhibited  prior  to  the  arrival  of  the  wave  at  the  S-IB  interface  implies  that  they 
are  not  a  consequence  of  interaction  with  the  specimen,  but  instead  are  imposed  upon  the 
specimen. 

Finally,  we  note  that  the  axial  velocity  and  axial  acceleration  of  the  specimen  at  the  S-IB 
interface  will  coincide  with  the  axial  velocity  and  acceleration  of  the  incident  bar  as  long  as  the 
specimen  and  bar  are  in  contact,  but  when  a  gap  forms  at  this  interface  the  velocities  and 
accelerations  of  the  bar  and  the  specimen  will  generally  be  different. 


6.  Simulations  with  the  Moony-Rivlin  Model  and  a  l-|ps  Rise  Time 


All  the  simulations  in  this  section  are  for  the  loading  wave  with  a  1-ps  initial  rise  time  and  the 
compressible  Mooney-Rivlin  model  for  the  specimen.  Conventions  pertaining  to  the  gaps  are 
discussed  in  section  2.5.  The  results  for  the  larger  radius  bar  (the  default  radius  Rb  =  12.8  mm  in 
this  report)  are  presented  in  subsections  6. 1-6.3.  These  results  are  compared  with  those  for  the 
smaller  radius  bar  in  subsection  6.4. 

6.1  Results  on  the  Centerline 

All  results  in  this  subsection  pertain  to  points  on  the  centerline  at  the  S-IB  or  the  S-TB  interface. 
We  begin  with  some  remarks  about  the  figures  and  a  few  general  observations,  and  then  proceed 
to  discuss  the  details  of  gap  formation  on  the  centerline. 


23  Of  course,  they  differ  quantitatively  by  a  factor  of  nearly  2  due  to  a  near  doubling  of  the  particle  velocity  on  reflection  of 
the  wave  at  the  S-IB  interface.  Also,  as  will  be  demonstrated  in  sections  6  and  7,  the  large  velocity  spikes  at  later  times  in  figures 
5  and  6  coincide  with  the  closing  of  gaps  on  the  centerline  at  the  S-IB  interface. 
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6.1.1  General  Remarks 

Figure  7  shows  the  time  histories  of  the  axial  displacement  of  the  incident  bar  (wIB)  and  the 
specimen  (us  m)  at  the  S-IB  interface,  as  well  as  the  history  of  the  corresponding  gap  size 
(As-ib  =  Ms- its  -  mib).  Similarly,  figure  8  shows  the  time  histories  of  the  axial  displacement  of  the 
transmission  bar  (mtb)  and  the  specimen  (ms-tb)  at  the  S-TB  interface,  as  well  as  the  history  of  the 
corresponding  gap  size  (AS-tb  =  mt b  -  ms-tb)-  In  all  of  the  figures  in  this  report,  the  displacements 
and  gaps  are  measured  in  microns.  The  histories  of  the  gaps  at  the  two  interfaces  are  compared 
in  figure  9.  Figure  10  enlarges  the  plot  in  figure  9  in  the  vicinity  of  the  time  at  which  the  gaps 
initially  form.  From  figures  7  to  10,  it  is  clear  that  gaps  first  form  at  both  interfaces  at  about 
165  ps.  These  gaps  open  and  close  several  times,  and  they  do  so  more  frequently  at  the  S-TB 
interface. 


Figure  7.  The  histories  of  the  axial  displacement  of  the  incident  bar  and  the  specimen  on  the  centerline  at 
the  S-IB  interface.  Also  shown  is  the  history  of  the  corresponding  gap  size  (blue). 
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Figure  8.  The  histories  of  the  axial  displacement  of  the  transmission  bar  and  the  specimen  on  the  centerline 
at  the  S-TB  interface.  Also  shown  is  the  history  of  the  corresponding  gap  size  (blue). 
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Figure  9.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  and  S-TB  interfaces. 
Also  shown  is  the  history  of  mean  axial  stretch  in  the  specimen  (blue). 
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Figure  10.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  and  S-TB  interfaces 
in  the  vicinity  of  the  time  at  which  they  first  form.  Also  shown  is  the  history  of  mean  axial 
stretch  in  the  specimen  (blue). 

The  largest  gaps  observed  were  43  pm  at  the  S-IB  interface  and  38  pm  at  the  S-TB  interface. 
Since  these  maximum  gap  sizes  are  smaller  than  the  initial  50-pm  axial  dimension  of  the 
specimen  elements,  the  issues  of  mesh  dependence  and  dependence  on  the  contact  algorithm 
naturally  arise.  These  issues  are  addressed  in  appendix  B.  The  results  described  there  support 
the  conclusion  that  the  gap  phenomenon  is  not  a  numerical  artifact  introduced  by  the  contact 
algorithm  or  stemming  from  an  insufficiently  fine  mesh.  Additional  support  for  the  latter 
conclusion  will  be  given  in  section  8  (linear  elastic  model  for  the  specimen),  where  gap  sizes  are 
found  to  exceed  the  mesh  size. 

Since  the  loading  wave  in  the  incident  bar  arrives  at  the  S-IB  interface  at  (or  slightly  before) 
t  =  140  ps,  whereas  the  gaps  first  form  about  25  ps  later,  gap  formation  does  not  coincide  with 
arrival  of  the  loading  wave.  Nor  does  gap  formation  coincide  with  the  first  wave  reflections 
from  either  interface.  Indeed,  since  the  travel  time  of  an  axial  longitudinal  wave  across  the 
specimen  is  about  0.7  ps  (see  table  2),  there  is  time  for  about  18  round  trips  of  this  wave  before 


29 


the  gaps  form.  In  any  case,  an  initially  compressive  longitudinal  wave  in  the  specimen  should 
reflect  from  the  transmission  bar  (and  subsequently  from  the  incident  bar)  as  a  compressive  wave 
due  to  the  higher  bar  impedance,  so  tensile  stresses  would  not  be  expected  to  develop  from  these 
wave  reflections.  On  the  other  hand,  it  is  clear  from  table  3  that  longitudinal  unloading  waves 
from  the  lateral  (stress-free)  surface  of  the  specimen  should  arrive  at  the  centerline  at  both 
interfaces  about  3-4  ps  after  the  loading  wave  in  the  bar  first  reaches  the  specimen,  but  this  does 
not  correlate  with  gap  formation  either.  Finally,  it  is  clear  from  table  2  that  25  ps  is  insufficient 
time  for  a  shear  wave  to  propagate  more  than  a  fraction  of  the  thickness  or  the  radius  of  the 
specimen;  consequently,  these  waves  cannot  influence  the  formation  of  gaps.  The  possible 
influence  of  wave  propagation  in  the  bars  on  the  formation  of  gaps  is  discussed  in  section  6.3. 

The  history  of  the  mean  axial  stretch  Az  (see  equation  8)  is  plotted  on  the  right  axis  in  figures 
9  and  10.  The  slope  of  the  axial  stretch  curve  in  figure  9  is  nearly  constant  for  t  >  200  ps,  and 
we  estimate  that  the  nominal  strain  rate  for  these  times  is  between  2470/s  and  2480/s.  This  is 
within  1%  of  the  desired  2500/s  strain  rate. 

Figure  1 1  plots  the  histories  of  the  axial  stress  azz  in  the  specimen  and  the  incident  bar  at  points 
near  the  centerline  and  the  specimen-bar  interface;  the  corresponding  gap  size  is  also  shown. 
Figure  12  is  an  analogous  plot  at  the  transmission  bar  interface.  The  stress  is  actually  measured 
at  the  centroids  of  the  specimen  and  bar  elements  adjacent  to  the  centerline  and  the  interface. 
Since  the  coordinates  of  these  centroids  do  not  coincide,  slight  differences  in  the  values  are  not 
unexpected.  Also,  recall  (section  2.6)  that  the  axial  stress  is  taken  to  be  negative  in  compression. 
When  a  gap  opens  between  the  faces  of  the  specimen  and  the  bar,  the  axial  stress  should  drop  to 
zero  on  both  faces,  regardless  of  the  size  of  the  gap.  The  plots  in  figures  1 1  and  12  are  consistent 
with  this  expectation,  keeping  in  mind  that  these  stresses  are  actually  measured  at  a  slight 
distance  away  from  the  faces. 
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Figure  1 1 .  Histories  of  the  axial  stress  (negative  in  compression)  and  the  gap  size  on  the  centerline  at 
the  S-IB  interface.  The  stress  is  measured  at  the  centroids  of  the  specimen  and  incident  bar 
elements  adjacent  to  the  centerline. 
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Figure  12.  Histories  of  the  axial  stress  (negative  in  compression)  and  the  gap  size  on  the  centerline  at 
the  S-TB  interface.  The  stress  is  measured  at  the  centroids  of  the  specimen  and 
transmission  bar  elements  adjacent  to  the  centerline. 

A  comparison  of  the  velocity  history  at  the  S-IB  interface  in  figure  5  with  the  displacement,  gap, 
stress  and  stretch  histories  in  figures  7-12  leads  to  the  following  additional  observations.  For 
convenience  of  discussion  we  have  divided  the  time  histories  into  several  stages. 

6.1.2  The  Initial  Loading  of  the  Specimen:  Stage  I  (140-165  ps) 

This  first  stage  covers  the  initial  loading  of  the  specimen  prior  to  gap  formation.  It  begins  with 
the  first  significant  increase  in  the  axial  velocity  of  the  incident  bar  at  the  S-IB  interface  at  about 
140  ps  (see  figure  5)24  and  ends  with  the  onset  of  the  gaps  at  165  ps.  From  the  axial  stretch 
histories  in  figures  9  and  10,  we  see  that  the  axial  stretch  decreases  from  1  to  about  0.97  over  this 
25-ps  time  interval;  hence,  the  nominal  axial  strain  increases  from  0%  to  about  3%. 

The  first  observable  increase  in  axial  stress  at  the  S-IB  interface  occurs  at  about  143  ps  (see 
figure  1 1).  The  stress  increases  in  magnitude  to  a  peak  value  of  -4.0  MPa  at  around  158  ps, 
which  is  7  ps  prior  to  the  formation  of  a  gap  at  165  ps.  After  this  peak,  the  stress  decreases  more 
rapidly  than  it  had  risen  and  drops  to  zero  at  about  164  ps.  Superimposed  on  this  behavior  are 

74 

"  However,  as  observed  at  the  end  of  section  5.3,  the  loading  wave  may  have  arrived  as  early  as  135  ps. 
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smaller  oscillations  with  a  period  of  1-2  ps;  similar  oscillations  can  be  seen  in  the  velocity 
history  at  the  S-IB  interface  (figure  5).  The  history  of  the  axial  stress  at  the  S-TB  interface  is 
qualitatively  and  quantitatively  similar,  with  the  times  delayed  by  about  1  ps  (see  figure  12).  In 
particular,  the  axial  stresses  at  the  two  interfaces  are  approximately  equal  throughout  this  stage, 
from  which  we  may  infer  that  ozz  is  approximately  uniform  through  the  thickness  of  the 
specimen.  In  the  experimental  literature,  this  condition  is  referred  to  as  dynamic  equilibrium.'25 

At  t  =  158  ps  (the  instant  of  the  stress  peak  in  figure  1 1),  the  axial  stretch  Az  is  about  0.99,  as  can 
be  inferred  from  figure  9.  This  stretch  corresponds  to  a  compressive  nominal  strain  of  0.01,  that 
is,  1%  axial  strain.  A  glance  at  the  stress-stretch  curve  for  uniaxial  stress  in  figure  A-l  reveals 
that  when  Az  =  0.99,  the  axial  stress  azz  is  a  very  small  fraction  of  1  MPa.  Indeed,  from  equation 
A-13  and  the  values  of  Ai  and  Ai  in  table  3,  we  find  that  ozz  «  -2.4  kPa  for  a  state  of  uniaxial 
stress  when  Az  =  0.99. 26  This  is  three  orders  of  magnitude  smaller  than  the  peak  value  of 
ozz  =  -4.0  MPa  in  figure  11,  which  (as  previously  noted)  occurs  at  Az  =  0.99.  Consequently,  the 
stress  state  at  the  instant  of  this  peak  (t  =  158  ps)  is  nowhere  close  to  the  desired  state  of  uniaxial 
stress.  In  fact,  as  established  in  the  next  paragraph,  at  this  instant  the  stress  state  on  the 
centerline  is  very  nearly  hydrostatic. 

Recall  from  section  2.6  that  at  points  on  the  centerline,  the  radial,  hoop  and  axial  directions  are 
principal  axes  of  stress  and  strain,  the  relations  14  hold  exactly,  and  equation  15  holds 
approximately  (for  small  volume  changes,  i.e.,  J  ~  1 ).  As  discussed  in  appendices  A-2  and  A-3, 
the  stress  difference  <jzz  -  o,T  on  the  centerline  is  given  exactly  by  equation  A- 14  when  7=1  and 
approximately  by  that  equation  when  J  ~  1 .  When  Az  =  0.99,  equation  A- 14  yields  the  value 
-2.4  kPa.27  Summarizing,  we  have  ozz  -  arr  ~  -2.4  kPa  and  dee  =  orr  exactly,  whereas 
gzz  =  -4.0  MPa.  Hence,  the  differences  between  the  principal  stresses  orr,  Gee,  and  gzz  are  either 
zero  or  at  least  three  orders  of  magnitude  smaller  than  the  stresses  themselves.  It  follows  that  all 
three  principal  stresses  are  nearly  equal  to  each  other  and  hence,  by  equation  11,  nearly  equal  to 
—p.  The  volumetric  strain  1-7  at  this  peak  pressure  is  approximately  0.001,  that  is,  0.1%. 28 

The  preceding  quantitative  estimates  apply  only  at  the  stress  peak  on  the  centerline.  However,  at 
neighboring  times  and  at  radial  locations  near  the  centerline  we  arrive  at  the  same  qualitative 
conclusion,  namely,  that  the  stress  state  is  nearly  hydrostatic  for  most  of  this  stage  (certainly 


25  It  should  be  kept  in  mind  that  the  stresses  in  figures  1 1  and  12  are  measured  on  the  centerline.  On  the  other  hand,  only  the 
mean  axial  stress  on  an  interface  can  be  inferred  from  experimental  data,  so  it  is  the  equality  of  the  mean  axial  stresses  at  the  two 
interfaces  that  is  implied  by  the  term  dynamic  equilibrium. 

25  As  discussed  in  appendix  A-3,  equation  A-13  holds  approximately  for  uniaxial  strain  when  /  « 1  (i.e.,  small  volume 
changes).  Also,  note  that  the  axial  strain  is  small  enough  that  the  exact  relation  in  equation  A-13  and  the  linear  elastic 
approximation  give  essentially  the  same  result. 

r)r7 

•i/  Since  the  right-hand  side  of  equation  A- 14  equals  the  right-hand  side  of  equation  A-13,  this  value  follows  from  the 
calculation  in  the  preceding  paragraph. 

28  This  follows  from  setting  p  equal  to  -a~  and  using  equation  A-21  with  K  =  4.0  GPa  (see  section  3.2),  which  yields 
1-7  =  4.0  MPa/4.0  GPa. 
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during  the  time  interval  147-163  jiis  in  which  la-1  exceeds  0.4  MPa),  and  the  axial  stress  exceeds 
that  for  a  uniaxial  stress  state  by  several  orders  of  magnitude. 

From  figure  5,  we  estimate  that  the  peak  slope  of  the  velocity-time  curve  during  Stage  I  also 
occurs  around  158  ps.29  Thus,  the  peak  in  the  axial  acceleration  of  the  incident  bar  at  the  S-IB 
interface  coincides  with  the  peak  in  the  axial  stress  on  the  centerline,  and  as  that  acceleration 
decreases  towards  zero  so  does  the  axial  stress.  The  axial  acceleration  of  the  incident  bar 
induces  an  axial  acceleration  in  the  specimen,  which  in  turn  induces  a  radial  acceleration  in  the 
specimen.  Therefore,  it  seems  likely  that  during  this  stage  the  value  of  the  axial  stress  and  the 
nearly  hydrostatic  stress  state  in  the  neighborhood  of  the  stress  peak  are  inertial  effects 
associated  with  the  accelerations  imposed  on  the  specimen  by  the  incident  bar.30  Furthermore, 
since  the  specimen  is  approximately  in  dynamic  equilibrium  (at  least  along  the  centerline),  these 
would  appear  to  be  primarily  radial  inertia  effects.  Song  et  al.  ( 6 )  also  observed  a  correlation 
between  axial  acceleration  and  axial  stress  in  the  specimen  at  early  times  in  their  SHPB  tests  and 
numerical  simulations  on  gel  specimens  and  likewise  concluded  that  this  was  a  radial  inertia 
effect.  This  conclusion  is  also  consistent  with  various  theoretical  analyses  of  inertial  effects 
applied  to  soft  materials  (26-33).  We  will  return  to  these  issues  in  section  7,  where  we  discuss 
the  case  with  a  25-ps  initial  rise  time. 

6.1.3  Gap  Opening  and  Closure  (Stages  II- VI) 

Stage  II  (165-280  ps):  This  stage  begins  with  the  opening  of  gaps  at  both  interfaces  and  ends 
with  the  first  closing  of  the  gap  at  the  S-IB  interface.  Over  this  1 15-ps  time  interval,  the  nominal 
axial  strain  increases  from  3%  to  32%  (see  figure  9). 

From  figure  10,  we  see  that  a  gap  begins  to  form  at  the  S-TB  interface  at  about  165.0  ps  and  at 
the  S-IB  interface  about  0.1  ps  later.  Prior  to  this  time,  As-ib  is  negative,  which  means  that  the 
specimen  and  incident  bar  have  interpenetrated.  The  maximum  amount  of  this  interpenetration  is 
0.06  pm,  which  is  0.12%  of  the  initial  axial  edge  length  of  a  specimen  element.  The  penalty- 
based  contact  algorithm  allowed  this  degree  of  interpenetration;  it  is  also  evident  in  several 
subsequent  figures.  This  interpenetration  begins  to  decrease  at  about  164.8  ps,  which  is  0.2  ps 
prior  to  the  instant  at  which  the  gap  forms  at  the  S-TB  interface.  The  axial  stretch  at  the  onset  of 
the  gaps  is  about  0.971.  The  corresponding  nominal  axial  strain  is  0.029,  that  is,  very  close  to 
3%. 

From  figure  5,  we  see  that  the  peak  value  of  the  axial  velocity  of  the  incident  bar  occurs  at  about 
164  ps.  The  axial  acceleration  of  the  bar  is  zero  at  this  instant  and  negative  for  some  time 
thereafter.  The  axial  stress  at  both  interfaces  also  reduces  to  zero  at  about  164  ps  (figures  11  and 
12).  Hence,  both  the  axial  stress  and  the  axial  acceleration  at  the  S-IB  interface  drop  to  zero  at 
about  the  same  time,  and  a  gap  forms  at  that  interface  about  1  ps  later,  that  is,  about  1  ps  after 


7Q 

Here,  we  disregard  the  smaller  oscillations  superimposed  on  the  main  curve. 
30  In  this  regard,  see  the  comments  in  the  next  to  last  paragraph  in  section  5.3. 
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the  incident  bar  begins  to  decelerate.  That  a  gap  forms  just  after  the  axial  stress  drops  to  zero  is 
not  surprising.  The  important  observation  here  is  that  this  occurs  when  the  axial  acceleration  of 
the  S-IB  interface  passes  through  zero.  We  will  return  to  this  point  in  section  7. 

With  the  exception  of  some  minor  oscillations,  the  incident  bar  continues  to  decelerate  until 
about  173  ps  and  then  re-accelerates  until  180  ps  (figure  5).  However,  the  S-IB  gap  does  not 
close  when  the  incident  bar  re-accelerates.  The  bar  velocity  goes  through  many  oscillations 
about  its  mean  value  of  about  2v0  before  the  S-IB  gap  first  closes  at  280  ps.  The  gap  at  the  S-TB 
interface  is  also  closed  at  280  ps,  although  it  has  opened  and  closed  several  times  prior  to  this. 
Nevertheless,  a  gap  exists  at  the  S-TB  interface  for  most  of  this  stage.  Observe  the  minor  spikes 
in  stress  at  the  S-TB  interface  when  the  gap  closes  there.  Note  that  the  magnitude  of  the  S-TB 
gap  never  exceeds  4  pm  in  this  stage,  whereas  the  S-IB  gap  reaches  a  maximum  of  38  pm  at 
about  226  ps. 

Stage  III  (280-286  ps):  The  gaps  at  both  interfaces  close  at  280  ps  and  remain  closed,  for  the 
most  part,  over  this  6-ps  time  interval.  The  S-IB  gap  does  re-open  momentarily  at  the  beginning 
of  this  stage,  and  sub-micron  gaps  re-open  momentarily  at  both  interfaces  near  the  end  of  the 
stage.  The  first  large  spike  in  the  velocity  at  the  S-IB  interface  (figure  5)  occurs  at  the  beginning 
of  this  stage,  and  large  spikes  in  the  stress  (figures  11  and  12)  occur  throughout  this  stage. 

Stage  IV  (286-324  ps):  The  gaps  re-open  and  remain  open  at  both  interfaces  over  this  30- ps 
time  interval.  The  axial  stresses  near  the  interfaces  again  drop  to  nearly  zero.  The  gaps  are  now 
roughly  the  same  size  at  both  interfaces.  The  S-IB  gap  reaches  a  peak  value  of  43  pm  at  308  ps; 
the  S-TB  gap  reaches  a  peak  value  of  38  pm  about  4  ps  later.  These  are  the  maximum  gap  sizes 
achieved  in  this  simulation.  During  this  stage,  the  nominal  axial  strain  increases  from  about  33% 
to  43%. 

Stage  V  (324-333  ps):  The  gaps  at  both  interfaces  close  and  remain  closed,  for  the  most  part, 
over  this  9-ps  time  interval.  This  gap  closure  results  in  the  second  large  spike  in  the  velocity  at 
the  S-IB  interface  (figure  5)  and  the  second  group  of  large  spikes  in  the  stress  (figures  11  and  12) 
as  the  specimen  “slaps”  the  bars. 

Stage  VI  (333-340  ps):  The  gaps  re-open  and  remain  open  at  both  interfaces.  The  axial  stresses 
near  the  interface  again  drop  to  nearly  zero.  At  340  ps,  the  axial  stretch  is  0.535  and  the  nominal 
axial  strain  is  46.5%.  Thus,  gaps  persist  to  large  strains  at  both  interfaces. 

6.2  Results  at  Other  Radial  Locations 

Figure  13  shows  deformed  mesh  plots  in  the  vicinity  of  the  specimen  at  236  ps  (Stage  II).  The 
centerline  lies  along  the  bottom  of  each  plot.  Thus,  figure  13a  shows  the  entire  (deformed) 
specimen  mesh  in  the  2-D  axisymmetric  simulation;  note  the  bulging  at  the  free  surface  of  the 
specimen  near  the  S-IB  interface.  A  close  examination  of  figure  13a  reveals  that  the  gap  at  the 
S-IB  interface  extends  outward  from  the  centerline  for  about  %  of  the  specimen  radius,  and  that  a 
gap  does  not  form  at  any  radial  location  along  the  S-TB  interface  at  this  instant.  The  S-IB  gap  is 
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clearly  visible  in  figure  13b  and  is  seen  to  be  approximately  uniform  in  size  for  much  of  its 
length.31  From  figure  9,  we  conclude  that  the  gap  at  the  S-IB  interface  is  about  37  pm  at  the 
centerline;  this  is  roughly  the  axial  dimension  of  the  deformed  specimen  elements.  The  axial 
stretch  at  this  instant  is  about  0.79,  corresponding  to  a  nominal  axial  strain  of  21%.  Figures 
8  and  9  also  confirm  that  there  is  no  gap  on  the  centerline  at  the  S-TB  interface  at  this  time. 


Figure  13.  Deformed  mesh  plots  in  the  vicinity  of  the  specimen  at  t  =  236  ps  and  an  axial  strain  of 

21%:  (a)  full  specimen  mesh;  (b)  enlargement  of  the  mesh  near  the  centerline,  where  the  gap  size 
is  about  37  pm. 

Figures  14  and  15  plot  the  gap  size  at  the  S-IB  interface  vs.  the  deformed  specimen  radius  r  at 
selected  times.  Figure  16  gives  the  same  information  at  the  S-TB  interface. 


31  The  latter  observation  is  consistent  with  the  results  in  figure  14. 
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Figure  14.  Gap  size  vs.  deformed  specimen  radius  at  the  S-IB  interface  at  selected  times. 


Figure  15.  Gap  size  vs.  deformed  specimen  radius  at  the  S-IB  interface  at  selected  times. 
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Figure  16.  Gap  size  vs.  deformed  specimen  radius  at  the  S-TB  interface  at  selected  times. 

The  selected  times  in  figure  14  span  Stage  II  and  also  include  the  instant  t  =  160  ps,  which  is 
prior  to  any  gap  formation.  Note  that  at  each  instant  after  160  ps,  the  largest  gap  occurs  at  a  non¬ 
zero  radial  coordinate,  that  is,  off  the  centerline.  The  largest  gap  size  during  this  stage  is  nearly 
42  pm,  which  occurs  at  230  ps  and  2  mm  from  the  centerline.  From  170  to  230  ps,  the  size  of 
the  gap  increases  with  time  and  the  gap  extends  along  most  of  the  specimen  radius,  that  is,  across 
most  of  the  face  of  the  specimen.  Over  this  time  interval,  the  nominal  axial  strain  increases  from 
4%  to  19%.  From  230  to  280  ps,  both  the  radial  extent  of  the  gap  and  the  maximum  size  of  the 
gap  decrease  with  time.  Over  this  time  interval,  the  nominal  axial  strain  increases  from  19%  to 
32%.  We  estimate  that  the  gap  extends  over  roughly  half  the  radius  of  the  specimen  at  t  =  265  ps 
and  a  nominal  axial  strain  of  28%. 

The  selected  times  in  figures  15  and  16  span  Stages  III- VI.  From  these  figures,  we  conclude  that 
the  radial  locations  where  the  S-IB  and  S-TB  gaps  in  Stage  IV  exceed  a  few  microns  are 
confined  to  the  vicinity  of  the  centerline.  The  largest  gap  sizes  during  this  time  interval  occur  at 
310  ps  and  are  nearly  equal:  43  pm  at  the  S-IB  interface  and  41  pm  at  the  S-TB  interface. 

6.3  Pressure  Contours 

Figures  17-23  display  pressure  contours  in  the  vicinity  of  the  specimen  at  specific  times  between 
152  and  166.4  ps.  Also  indicated  at  the  top  of  each  contour  plot  is  the  mean  axial  stretch  Xz  at 
the  selected  time.  Each  contour  plot  shows  only  a  small  portion  of  the  length  of  the  bars  but  does 
span  the  bars  from  the  centerline  (bottom)  to  the  stress-free  lateral  surface  (top).  The  incident 
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bar  is  on  the  left  and  the  transmission  bar  is  on  the  right.  Recall  that  the  radius  of  the  bars  is 
slightly  more  than  twice  the  initial  radius  of  the  specimen.  Hence,  the  large,  rectangular,  white 
region  between  the  bars  in  these  contour  plots  is  just  empty  space  above  the  specimen.  The 
specimen-bar  interfaces  are  just  barely  distinguishable  as  the  thin,  light  lines  below  this  region. 
The  length  scale  on  these  figures  can  be  inferred  from  the  thickness  of  the  specimen,  which  is 
initially  1.45  mm  and  does  not  change  substantially  over  the  time  interval  considered  here. 

Since  there  is  no  load  on  the  face  of  a  bar  where  it  overhangs  the  specimen,  the  stress  on  that  part 
of  the  surface  is  zero,  and  since  the  bar  faces  remain  very  nearly  planar,  we  may  conclude  that 
ozz,  <Jrz,  and  Go-  are  essentially  zero  on  the  overhanging  faces.  But  the  stress  tensor  a  need  not  be 
identically  zero  at  points  on  the  overhanging  surface,  so  the  pressure  need  not  be  zero  there. 
Indeed,  the  radial  and  hoop  stresses  in  the  bars  could  be  nonzero  at  some  points  on  the 
overhanging  surface,  which  would  result  in  a  nonzero  pressure  there;  cf.  equation  ll.32 
Similarly,  at  those  places  and  times  at  which  the  specimen  and  the  bars  are  in  contact,  the 
pressure  need  not  be  continuous  across  the  interface.  Since  the  contact  is  frictionless,  the  only 
nonzero  component  of  the  traction  vector  at  the  specimen-bar  interface  is  the  axial  stress  ozz,  so  it 
is  this  stress  component  that  must  be  continuous  across  the  interface  when  there  is  contact. 
Likewise,  it  is  the  axial  stress  (but  not  necessarily  the  pressure)  that  must  be  zero  on  the 
specimen  and  bar  faces  at  points  where  a  gap  forms. 

Recall  that  pressure  is  taken  to  be  positive  in  compression  (equation  11).  The  light  green  areas  in 
the  contour  plots  are  regions  of  positive  pressure,  whereas  the  dark  green  areas  are  regions  of 
negative  pressure.  The  dark  (and  typically  jagged)  curves  bounding  these  two  green  regions  are 
the  zero-pressure  level  curves.  Clearly,  the  pressure  field  in  the  specimen  and  the  bars  is  an 
extremely  complex  function  of  position  and  time.  Reflection  of  waves  from  the  interfaces  and 
the  free  surfaces  and  oscillations  induced  in  the  specimen  and  bar  by  the  passage  of  these  waves 
make  it  difficult  to  sort  out  cause  and  effect  here,  particularly  in  view  of  the  fact  that  the  pressure 
is  only  part  of  the  full  stress  state. 

Figures  17-20  display  pressure  contours  at  specific  times  between  152  and  164.8  ps,  which  is 
prior  to  the  formation  of  a  gap.  Negative  (tensile)  pressures  develop  in  the  bars  prior  to 
developing  in  the  specimen.  Even  at  the  early  time  of  152  ps  (figure  17),  there  are  small  regions 
of  negative  pressure  in  portions  of  the  bars  that  overhang  the  specimen.  These  tensile  regions  are 
most  likely  generated  by  unloading  waves  from  the  stress-free  portion  of  the  bar  faces.  By 
163.4  ps  (figure  18),  negative  pressure  has  developed  in  the  incident  bar  in  a  region  bordering 
the  centerline,  about  6-8  specimen  lengths  back  from  the  S-IB  interface.  From  164.2  tol64.6  ps 
(figures  19  and  20),  the  pressure  in  the  incident  bar  oscillates  (with  axial  distance  from  the 
specimen)  between  positive  and  negative  values. 


32  Likewise,  since  the  low  amplitude  waves  in  the  bars  should  induce  little  bulging  of  the  stress-free  lateral  surface  of  the 
bars,  a„,  a9r,  and  azr  are  essentially  zero  on  this  surface,  although  the  pressure  need  not  be  zero  there. 
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Figure  17.  Pressure  contours  in  the  vicinity  of  the  specimen  at  selected  times. 
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Figure  18.  Pressure  contours  in  the  vicinity  of  the  specimen  at  selected  times. 
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Negative  pressure  does  not  develop  in  the  specimen  until  164.4  ps  (figure  20),  although  it  is 
limited  to  very  small  regions.  By  164.8  ps,  most  of  the  specimen  and  the  regions  of  the  bars 
adjacent  to  it  are  in  a  state  of  negative  pressure.  Gap  formation  begins  after  164.8  ps  (see  section 
6.1.3).  From  165.0  ps  onward,  (figures  21-23),  the  regions  of  the  bars  adjacent  to  the  specimen 
are  primarily  in  a  state  of  negative  pressure.  However,  the  sign  of  the  pressure  in  the  specimen 
oscillates  (with  time)  from  primarily  negative  (164.8-165.2  ps)  to  primarily  positive 
(165.4-165.8  ps)  and  back  again  to  primarily  negative  (166.0-166.4  ps).  In  view  of  this 
oscillation  in  the  sign  of  the  pressure,  it  is  possible  that  the  axial  stress  also  oscillates  between 
tensile  and  compressive  values  at  points  in  the  interior  of  the  specimen,33  although  this  cannot  be 
confirmed  from  the  pressure  contours. 


33  Recall  that  tensile  axial  stresses  cannot  be  supported  at  the  specimen  and  bar  faces. 
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Figure  23.  Pressure  contours  in  the  vicinity  of  the  specimen  at  selected  times. 

Finally,  note  that  at  152  and  162  jus  (figure  17),  the  pressure  on  the  centerline  in  (and  near)  the 
specimen  is  around  1  MPa.  This  is  consistent  with  the  conclusion  that  the  stress  state  on  the 
centerline  at  these  instants  is  nearly  hydrostatic  (see  section  6.1.2).  Indeed,  from  figures 
1 1  and  12,  we  see  that  azz  is  about  -1  MPa  at  these  times. 

6.4  Effects  of  Pressure  Bar  Radius  on  Gap  Formation 

The  results  in  section  6.3  suggest  that  the  large  overhang  of  the  bars  (particularly,  the  incident 
bar)  may  contribute  to  the  formation  of  gaps  in  the  specimen.  In  order  to  explore  this  hypothesis, 
we  reduced  the  bar  radius  i?B  from  12.8  to  9.5  mm  (see  section  2.1).  The  histories  of  the  gap  size 
on  the  centerline  for  both  cases  are  compared  in  figure  24.  It  is  clear  that  the  reduction  in  the  bar 
radius  substantially  reduced  the  size  of  the  gaps  at  both  interfaces  at  later  times. 
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Figure  24.  Histories  of  the  gap  size  on  the  centerline  at  the  S-IB  interface  (top)  and  S-TB  interface 
(bottom)  for  two  different  bar  radii. 


48 


Somewhat  surprisingly,  however,  figure  24  reveals  that  this  reduction  in  the  bar  radius  increased 
both  the  size  and  the  duration  of  the  initial  gaps  at  both  interfaces.  Also,  for  the  smaller  diameter 
bar,  the  gaps  at  both  interfaces  form  about  2  ps  earlier.  These  observations  suggest  a  stronger 
interaction  with  the  stress-free  lateral  surfaces  of  the  bars,  which  are  now  closer  to  the  specimen. 

The  effect  of  the  bar  radius  on  the  formation  of  gaps  for  annular  specimens  is  discussed  in 
appendix  C-2,  but  otherwise  all  results  that  follow  are  for  the  larger  radius  Rb  =  12.8  mm. 


7.  Simulations  with  the  Moony-Rivlin  Model  and  a  25-|iis  Rise  Time 


All  the  simulations  in  this  section  are  for  the  loading  wave  with  a  25 -ps  initial  rise  time  and  the 
compressible  Mooney-Rivlin  model  for  the  specimen.  For  this  case,  the  axial  velocity  histories 
at  four  locations  on  the  centerline  of  the  incident  bar  are  given  in  figures  4  and  6.  As  observed  in 
section  5.3,  the  velocity  history  at  the  S-IB  interface  for  this  case  (figure  6)  is  substantially 
smoother  than  the  corresponding  velocity  history  for  the  1-ps  initial  rise  time  (figure  5).  Recall 
that  the  final  rise  time  for  the  1-ps  case  is  24  ps.  For  the  25-ps  initial  rise  time,  the  final  rise 
time  has  increased  to  about  41  ps.  The  purpose  of  this  section  is  to  examine  the  effects  of  this 
“pulse  shaping”  on  the  formation  of  gaps. 

7.1  Results  on  the  Centerline 

All  results  in  this  subsection  pertain  to  points  on  the  centerline  at  the  S-IB  or  S-TB  interface.  We 
begin  with  some  remarks  about  the  figures  and  a  few  observations,  and  then  proceed  to  discuss 
the  details  of  gap  formation  on  the  centerline. 

7.1.1  General  Remarks 

Figure  25  compares  the  gap  size  histories  on  the  centerline  at  the  S-IB  interface  for  the  1-  and 
25-ps  initial  rise  times;  figure  26  does  the  same  at  the  S-TB  interface.  Clearly,  the  pulse  shaping 
has  substantially  reduced  the  size  and  the  duration  of  the  gaps  at  both  interfaces. 
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Figure  25.  Histories  of  the  gap  size  on  the  centerline  at  the  S-IB  interface  for  two  loading  waves  with 
different  initial  rise  times. 
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Figure  26.  Histories  of  the  gap  size  on  the  centerline  at  the  S-TB  interface  for  two  loading  waves  with 
different  initial  rise  times. 

The  gap  size  histories  for  the  25 -ps  initial  rise  time  are  more  easily  discernable  as  the  blue  curves 
in  figure  27  (S-IB  interface)  and  figure  28  (S-TB  interface).  The  largest  gap  at  the  S-IB  interface 
is  now  only  2.45  pm,  while  the  largest  gap  at  the  S-TB  interface  is  now  only  1.27  pm.  The 
sensitivity  of  the  gap  size  histories  to  the  mesh  size  and  the  contact  algorithm  parameters  is 
examined  in  appendix  B.  The  results  described  there  support  the  conclusion  that  the  gap 
phenomenon  is  not  a  numerical  artifact  introduced  by  the  contact  algorithm  or  stemming  from  an 
insufficiently  fine  mesh. 
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Figure  27.  Histories  of  the  axial  stress  (negative  in  compression)  and  the  gap  size  on  the  centerline  at  the  S-IB 
interface.  The  stress  is  measured  at  the  centroids  of  the  specimen  and  incident  bar  elements 
adjacent  to  the  centerline. 
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Figure  28.  Histories  of  the  axial  stress  (negative  in  compression)  and  the  gap  size  on  the  centerline  at  the 
S-TB  interface.  The  stress  is  measured  at  the  centroids  of  the  specimen  and  incident  bar 
elements  adjacent  to  the  centerline. 

Figure  27  also  plots  the  histories  of  the  axial  stress  ozz  in  the  specimen  and  the  incident  bar  at 
points  near  the  centerline  and  the  specimen-bar  interface.  Figure  28  contains  an  analogous  plot 
at  the  transmission  bar  interface.  Recall  that  the  stress  is  actually  measured  at  the  centroids  of 
the  specimen  and  bar  elements  adjacent  to  the  centerline  and  the  interface.  Since  the  coordinates 
of  these  centroids  do  not  coincide,  slight  differences  in  the  values  are  possible.  Also,  recall  that 
the  axial  stress  is  taken  to  be  negative  in  compression.  As  expected,  the  axial  stress  is  zero  or 
nearly  zero34  when  a  gap  exists,  even  when  the  gap  is  submicron  in  size.  Sharp  stress  spikes 
occur  when  the  gap  closes  and  the  specimen  “slaps”  the  bar. 

We  did  not  plot  the  mean  axial  stretch  history  for  the  case  considered  in  this  section.  All  stretch 
(and  corresponding  strain)  values  cited  below  were  estimated  from  the  axial  stretch  history  for  a 


34  Since  the  element  centroids  do  not  lie  on  the  interface,  the  axial  stress  may  have  a  small  but  nonzero  value  at  the  centroid 
even  if  it  is  zero  on  the  face  of  the  element  due  to  the  formation  of  a  gap. 
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linear  elastic  specimen  with  a  25-ps  initial  rise  time;35  see  figure  36  in  section  8.2.  The  slope  of 
the  axial  stretch  curve  in  that  figure  is  constant  for  t  >  190  ps,  and  we  estimate  that  the  nominal 
strain  rate  for  these  times  is  about  2490/s.36  This  is  within  0.5%  of  the  desired  2500/s  strain  rate. 

A  comparison  of  the  velocity  history  at  the  S-IB  interface  in  figure  6  with  the  gap  and  stress 
histories  in  figures  27  and  28  leads  to  the  following  additional  observations.  For  convenience  of 
discussion  we  have  divided  the  time  histories  into  four  stages.  The  first  two  stages  are  analogous 
to  those  for  the  1-ps  initial  rise  time  (see  sections  6.1.2  and  6.1.3). 

7.1.2  The  Initial  Loading  of  the  Specimen:  Stage  I  (143-184  ps) 

This  first  stage  covers  the  initial  loading  of  the  specimen  prior  to  gap  formation.  It  begins  with 
the  first  significant  increase  in  the  axial  velocity  of  the  incident  bar  at  the  S-IB  interface  at  about 
143  ps  (see  figure  6)  and  ends  with  the  onset  of  the  first  gap  at  184  ps.  The  axial  stretch 
decreases  from  1  to  about  0.95  over  this  time  41-ps  interval;  hence,  the  nominal  axial  strain 
increases  from  0%  to  about  5%. 

The  first  observable  increase  in  axial  stress  occurs  at  about  149  ps  at  the  S-IB  interface  and 
150  ps  at  the  S-TB  interface.  The  stress  in  the  specimen  increases  to  a  sharp  peak  of 
-1.6  and  -1.5  MPa  at  the  S-IB  and  S-TB  interfaces,  respectively,  at  about  164  ps.  Thereafter, 
the  stress  at  both  interfaces  decreases  to  a  broad  local  minimum  centered  at  about  173  ps,  then 
increases  again  to  a  local  maximum  around  177-178  ps,  and  finally  decreases  to  zero  at  184  ps. 
The  axial  stress  at  this  second  (minor)  peak  is  -1.1  and  -1.0  MPa  at  the  S-IB  and  S-TB 
interfaces,  respectively.  Smaller  oscillations  are  superimposed  on  this  general  behavior. 

The  axial  stretch  at  164  ps,  the  time  of  the  first  (major)  peak  in  stress,  is  about  0.99.  Since  this  is 
the  same  value  of  the  stretch  at  the  stress  peak  for  the  1-ps  initial  rise  time  (section  6.1.2),  we 
conclude  that  at  this  instant,  ozz  -  arr  «-2.4  kPa  and  aee  =  a,r.  However,  as  noted  above, 
ozz  w  -1.5  MPa.  Hence,  the  differences  between  the  principal  stresses  <jrr,  aee,  and  azz  are  either 
zero  or  nearly  three  orders  of  magnitude  smaller  than  the  stresses  themselves.  It  follows  that  all 
three  principal  stresses  are  nearly  equal  to  each  other  and,  hence,  to  —p.  That  is,  the  stress  state  is 
nearly  hydrostatic. 

While  the  preceding  quantitative  estimates  apply  only  at  the  stress  peak,  similar  calculations  at 
neighboring  times  lead  to  the  same  qualitative  conclusion:  the  stress  state  is  nearly  hydrostatic 
for  most  of  this  stage,  as  opposed  to  the  desired  state  of  uniaxial  stress. 

A  close  examination  of  the  axial  velocity  history  of  the  incident  bar  at  the  specimen  interface 
(figure  6)  reveals  that  the  extrema  in  the  axial  stress  in  the  specimen  at  164,  173,  and  177-178  ps 
coincide  (at  least  approximately)  with  inflection  points  in  the  velocity  vs.  time  curve,  that  is,  with 

35  Recall  that  the  mean  axial  stretch  is  computed  from  equation  8.  Because  of  the  low  impedance  of  the  specimen  for  both 
models,  the  displacements  and  hence  the  mean  stretches  for  these  two  cases  should  be  close. 

36  get  essentially  the  same  value  using  the  approximation  for  the  nominal  strain  rate  in  equation  10,  together  with 
Vjb  =  1.99  Vq,  which  is  the  mean  value  for  the  velocity  of  the  incident  bar  at  later  times,  estimated  from  figure  6. 
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extrema  in  the  axial  acceleration  vs.  time  curve.  Similarly,  as  the  acceleration  decreases  towards 
zero  after  178  ps,  so  does  the  axial  stress.  Thus,  just  as  for  the  case  of  the  1-jlis  initial  rise  time 
(section  6.1.2),  it  seems  likely  that  during  this  stage,  the  value  of  the  axial  stress  and  the  nearly 
hydrostatic  stress  state  in  the  neighborhood  of  the  stress  peak  are  inertial  effects  associated  with 
the  accelerations  imposed  on  the  specimen  by  the  incident  bar.37  Furthermore,  since  the  axial 
stresses  at  either  face  are  nearly  equal,  we  may  conclude  that  the  axial  stress  is  approximately 
uniform  through  the  thickness  of  the  specimen  (at  least  on  the  centerline),  so  these  would  appear 
to  be  primarily  radial  inertia  effects. 

The  preceding  conclusion  is  reinforced  by  comparing  the  peak  stresses  and  accelerations  for  the 
two  different  rise  times.  The  peak  axial  stress  at  the  S-IB  interface  for  the  1-ps  case  (-4.0  MPa) 
is  about  2Vi  times  the  peak  stress  (-1.6  MPa)  for  the  25-ps  case.  From  the  velocity  histories  in 
figures  5  and  6,  we  can  obtain  rough  estimates  for  the  axial  accelerations  of  the  incident  bar  at 
the  time  of  these  peak  stresses.  We  find  that  the  acceleration  for  the  1-ps  case  is  between  2Vi  to 
3  times  that  for  the  25-ps  case.  Thus  the  peak  axial  stress  in  this  early  stage  increases  in  (rough) 
proportion  to  the  peak  axial  acceleration.  Song  et  al.  ( 6 )  reached  a  similar  conclusion  from  their 
numerical  simulations  of  SHPB  tests  on  gel  rubber  specimens.  This  conclusion  is  also  consistent 
with  various  theoretical  analyses  of  inertial  effects  for  soft  materials  (26-33). 

7.1.3  Gap  Opening  and  Closure  (Stages  II-IV) 

Stage  II  (184-218  ps):  This  stage  begins  with  the  opening  of  the  first  gap,  which  occurs  at  the 
S-TB  interface,  and  ends  with  the  first  closing  of  the  gap  at  the  S-IB  interface.  During  this  34- ps 
time  interval,  the  axial  strain  increases  from  5%  to  about  14%.  Note  that  the  duration  of  this 
stage  is  about  a  third  of  the  duration  of  Stage  II  for  the  case  of  a  1  -ps  initial  rise  time,  and  gap 
formation  occurs  about  20  ps  later  than  it  did  for  the  1-ps  case. 

The  axial  stress  at  both  interfaces  drops  to  zero  at  about  184  ps.  The  gap  at  the  S-TB  interface 
forms  at  this  time,  whereas  the  gap  at  the  S-IB  interface  forms  about  2  ps  later.  Analogous  to  the 
case  for  the  1-ps  initial  rise  time,  the  contact  algorithm  allowed  a  0.10- pm  interpenetration  at  the 
S-IB  interface  prior  to  gap  formation.  This  interpenetration  begins  to  decrease  at  about  184  ps. 
From  figure  6,  we  see  that  the  peak  in  the  axial  velocity  of  the  incident  bar  at  the  S-IB  interface 
also  occurs  at  184  ps.  The  axial  acceleration  of  the  bar  is  zero  at  this  instant  and  then  becomes 
negative  for  some  time  thereafter. 

This  stage  ends  with  the  first  closing  of  the  gap  at  the  S-IB  interface  at  about  218  ps.  The  gap  at 
the  S-TB  interface  is  also  closed  at  this  instant,  although  it  has  opened  and  closed  several  times 
previously.  Nevertheless,  it  is  open  for  most  of  this  stage.  The  largest  gap  size  at  the  S-IB 
interface  (but  not  at  the  S-TB  interface)  is  attained  during  this  stage:  about  2.45  pm  at  200  ps. 


37  In  this  regard,  see  the  comments  in  the  next  to  last  paragraph  in  section  5.3. 
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Stage  III  (218-302  jos):  This  stage  begins  with  the  first  closing  of  the  gap  at  the  S-IB  interface, 
which  produces  large  spikes  in  the  axial  stress  at  both  interfaces  (figures  27  and  28)  and  a  large 
spike  in  the  axial  velocity  at  the  S-IB  interface  (figure  6).  Subsequently,  the  gaps  at  both 
interfaces  open  and  close  many  times,  again  with  large  spikes  in  stress  and  velocity 
corresponding  to  gap  closure.  This  stage  ends  with  the  permanent  closing  of  the  gaps  at  both 
interfaces  at  about  302  ps.  During  this  84  ps  time  interval,  the  nominal  axial  strain  increases 
from  14%  to  34%.  The  largest  gap  at  the  S-TB  interface  occurs  during  this  stage:  1.27  pm  at 
t  =  247  ps. 

Stage  IV  (302-340  ps):  The  gaps  at  both  interfaces  close  and  remain  closed  during  this  final 
stage.  The  axial  stress  oscillates  (without  spikes)  about  a  steadily  increasing  mean  value.  There 
is  no  analog  of  this  stage  for  the  case  of  a  1-ps  initial  rise  time,  at  least  for  times  less  than  340  ps. 

During  Stage  IV,  the  axial  stretch  Az  decreases  monotonically  from  0.66  to  0.57,  and  hence  the 
nominal  axial  strain  increases  from  34%  to  43%.  From  equation  A- 3  and  the  values  of  A  i  and  A  2 
in  table  3,  we  find  that  for  this  range  of  Az,  -a,,  would  increase  monotonically  from  about 
1 17  kPa  to  175  kPa  if  the  stress  state  were  uniaxial  (see  also  figure  A-l).  However,  from  figures 
27  and  28  we  see  that  in  fact  -ozz  oscillates  between  0  and  1.1  MPa  during  this  stage,  with  an 
average  value  around  0.5  MPa  =  500  kPa.  Hence,  the  specimen  is  not  in  a  state  of  uniaxial  stress 
during  this  stage.  Indeed,  using  equation  A- 14  and  arguing  as  in  section  6.1.2,  we  conclude  that 
for  most  Stage  IV,  the  stress  state  is  closer  to  hydrostatic  than  it  is  to  uniaxial. 

7.2  Gap  Size  vs.  Radius  at  Selected  Times 

Figure  29  plots  the  gap  size  vs.  deformed  specimen  radius  at  the  S-IB  interface  at  selected  times. 
There  is  no  gap  at  any  radial  location  for  times  less  than  or  equal  to  185  ps.  Also,  note  the 
0.10  pm  interpenetration  allowed  by  the  contact  algorithm  at  the  S-IB  interface  prior  to  gap 
formation. 
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Figure  29.  Gap  size  vs.  deformed  specimen  radius  at  the  S-IB  interface  at  selected  times. 

At  190  (as,  the  gap  extends  along  most  of  the  face  of  the  specimen.  Thereafter,  the  radial  extent 
of  the  gap  decreases  with  time.  By  215  ps,  the  gap  extends  only  3  mm  out  from  the  centerline, 
which  is  slightly  less  than  half  the  radius.  The  maximum  gap  size  increases  with  time  up  to 
200  ps  and  then  decreases  with  time  until  about  220  ps,  when  the  gap  appears  to  close 
temporarily.  Also,  note  that  at  each  instant  from  190-215  ps,  the  maximum  size  of  the  gap 
occurs  off  the  centerline.  From  230  ps  onward,  gap  opening  and  closing  exhibits  a  more 
complicated  radial  and  temporal  behavior.  The  results  in  figure  29  are  consistent  with  the  those 
on  the  centerline  in  figures  25  and  27. 

7.3  Force  History  in  the  Transmission  Bar 

Figure  30  is  somewhat  similar  to  figure  28;  in  particular,  the  history  of  the  gap  size  on  the 
centerline  at  the  S-TB  interface  is  repeated  in  figure  30.  But  instead  of  plotting  the  axial  stress  in 
the  transmission  bar  at  a  point  very  near  the  S-TB  interface,  figure  30  plots  the  total  axial  force 
Fz  acting  on  the  cross-section  of  the  transmission  bar  at  an  axial  location  2.52  mm  from  the 
specimen-bar  interface.  This  is  about  25  bar  elements  back  from  the  interface.  Also,  note  that  Fz 
has  been  taken  positive  in  compression,  whereas  azz  in  figure  28  is  positive  in  tension.38 


38  When  comparing  figures  28  and  30,  keep  in  mind  that  the  colors  of  the  stress-force  curves  and  the  gap  size  curves  have 
been  switched,  as  have  the  left-right  locations  of  the  vertical  axes. 
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Figure  30.  History  of  the  total  force  (positive  in  compression)  at  the  stress  gage  location  in  the 

transmission  bar  (2.52  mm  from  the  specimen  interface)  compared  with  the  history  of  the 
gap  size  on  the  centerline  at  the  S-TB  interface. 

On  taking  the  difference  in  the  sign  conventions  into  account,  we  see  (not  unexpectedly)  that  the 
force  history  in  figure  30  is  qualitatively  similar  to  the  stress  history  in  figure  28.  In  particular, 
the  force  is  positive  from  150  to  186  ps,  corresponding  to  the  Stage  I  loading  discussed  in 
section  7.1.2,  and  the  timing  of  the  local  maximum  and  minimum  values  of  Fz  and  azz  coincide 
during  this  stage.  A  gap  opens  and  closes  repeatedly  at  the  centerline  starting  at  184  ps; 
subsequently,  the  total  force  oscillates  with  repeated  dips  below  zero  (i.e.,  a  tensile  force).  This 
continues  until  about  260  ps,  after  which  the  force  remains  positive  and  oscillates  about  a 
steadily  increasing  mean  value. 

When  comparing  the  force  and  gap  histories  in  figure  30,  it  should  be  kept  in  mind  that  the  total 
force  is  the  cumulative  effect  of  the  stress  at  all  radial  locations  in  the  transmission  bar,  whereas 
the  gap  has  been  measured  only  on  the  centerline.  Nevertheless,  from  184  to  302  ps  (Stages  II 
and  III  in  section  7.1.3),  the  local  maximum  and  minimum  values  of  the  total  force  coincide,  for 
the  most  part,  with  the  closure  and  opening,  respectively,  of  a  gap  on  the  centerline.  In 
particular,  the  times  at  which  the  total  force  becomes  tensile  coincide  with  gaps,  for  the  most 
part.  It  seems  reasonable  to  conclude  that  when  a  gap  exists  on  the  centerline  at  the  S-TB 
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interface,  it  extends  over  a  significant  portion  of  the  specimen  face,  although  we  do  not  have  any 
plots  (analogous  to  figure  29)  to  verify  this. 

The  force  in  the  transmission  bar  in  figure  30  is  measured  at  an  axial  location  2.52  mm  from  the 
specimen-bar  interface.  This  location  was  selected  because  it  is  typical  for  the  placement  of  a 
quartz  stress  gage  in  the  bar  in  SHPB  tests  (Casern  et  al.  [25]).  Hence,  we  refer  to  this  axial 
location  as  the  gage  location,  even  though  a  stress  gage  was  not  included  in  the  simulation.  The 
signal  from  a  quartz  gage  is  proportional  to  the  total  force  applied  to  it,  regardless  of  whether  or 
not  the  axial  stress  in  the  transmission  bar  is  radially  uniform  (25).  With  appropriate 
calibrations,  the  force  measured  by  a  quartz  gage  should  also  coincide  with  the  force  inferred 
from  semi-conductor  strain  gage  measurements  midway  down  the  bar;  see  Moy  et  al.  (5)  for  the 
case  of  a  ballistic  gelatin  specimen. 

In  view  of  the  proximity  of  the  gage  location  to  the  S-TB  interface,  we  expect  that  the  force  on 
the  transmission  bar  at  the  specimen  interface  may  be  approximated  by  the  force  at  the  gage 
location,  with  some  qualifications  (see  below).  Hence,  keeping  in  mind  the  sign  conventions  and 
letting  A0  denote  the  original  cross-sectional  area  of  the  specimen,  we  see  that  -FJAo  provides 
an  estimate  for  the  average  or  mean  value  of  the  axial  component  of  the  nominal  stress  in  the 
specimen.39  Similarly,  dzz ,  the  mean  value  of  the  true  axial  stress  ozz  in  the  specimen,  is 

estimated  by  dividing  -F.  by  the  deformed  cross-sectional  area  of  the  specimen,  which  in  turn  is 
approximated  by  Ao/zU  ,40  so  that  azz  «  -XzFzIAq.  These  estimates  for  the  mean  value  of  the 

nominal  and  true  axial  stress  in  the  specimen  are  consistent  with  those  used  in  reporting  data  in 
the  SHPB  literature. 

Let  us  apply  this  estimate  for  dzz  to  the  first  peak  in  F-  at  164-165  ps,  which  is  prior  to  gap 

formation.  The  force  at  this  instant  is  25  N,  Az  is  about  0.99  (see  section  7.1.2),  and  the 
undeformed  specimen  radius  is  6.35  mm,  which  gives  azz  «-0.20  MPa  =  -200  kPa.  This  is 

almost  two  orders  of  magnitude  above  the  value  of  -2.4  kPa  that  would  correspond  to  a  state  of 
uniaxial  stress  at  this  axial  stretch.  On  the  other  hand,  the  axial  stress  azz  on  the  centerline  at  the 
S-TB  interface  at  this  instant  is  -1.5  MPa  (section  7.1.2),  which  is  7-8  times  the  estimated  mean 
stress.  Hence,  if  this  estimate  for  the  mean  stress  is  valid  (at  least  approximately),  then  the  axial 
stress  must  be  radially  non-uniform,  with  higher  values  near  the  centerline.  This  qualitative 
behavior  is  consistent  with  theoretical  analyses  of  radial  inertia  effects  for  soft  materials  (7,  26, 
27,  32,  33)  and  with  the  results  of  other  numerical  simulations  of  SHPB  tests  (6,  34). 

For  early  times,  the  force  history  in  figure  30  is  qualitatively  similar  to  the  mean  nominal  stress 
histories  on  solid  specimens  in  figure  3  of  Moy  et  al.  (5)  and  figure  5  of  Song  et  al.  (6).  The 
latter  two  figures  report  experimental  SHPB  data  on  a  ballistic  gelatin  specimen  and  a  gel  rubber 
specimen,  respectively.  All  three  figures  show  a  large  “inertial  spike”  at  small  strains  and  early 


39  This  is  also  referred  to  as  the  engineering  stress.  It  is  the  zz-component  of  the  1st  Piola-Kirchhoff  stress  tensor. 

This  approximate  relation  is  exact  when  the  axial  stretch  in  the  specimen  is  uniform  and  there  is  no  volume  change. 
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times  (165  ps  in  figure  30).  In  Moy  et  al.  (5)  and  Song  et  al.  ( 6 )  it  is  demonstrated  that  this 
inertial  spike  can  be  substantially  reduced  by  hollowing  out  the  specimen,  the  rational  being  that 
radial  inertia  contributions  to  the  axial  stress  are  largest  in  the  center  of  the  specimen  (cf.  the 
preceding  paragraph).  Unlike  figure  30,  the  experimental  data  in  Moy  et  al.  (5)  does  not  exhibit 
the  large  oscillations  or  any  evidence  of  gap  formation  after  the  inertial  spike.  This  is  not 
surprising  in  view  of  the  fact  that  the  pulse  shaping  in  those  tests  resulted  in  a  much  smoother 
loading  wave  than  the  one  used  in  our  simulation.  On  the  other  hand,  the  entire  force  history  in 
figure  30  is  qualitatively  similar  to  nominal  stress  history  for  the  solid  specimen  in  figure  5  of 
Song  et  al.  (<5),  for  which  no  pulse  shaping  was  used.  These  similarities  include  the  oscillations 
and  the  dips  to  negative  (i.e.,  tensile)  forces  following  the  inertial  spike.  Thus  the  data  in  Song  et 
al.  ( 6 )  may  be  exhibiting  indirect  evidence  of  gap  formation  at  the  S-TB  interface,  although  the 
possibility  of  gaps  was  not  discussed  in  that  paper.41 

Based  on  the  discussion  above  and  on  the  previous  discussion  of  the  histories  of  the  stress, 
velocity,  and  gap  size  on  the  centerline  for  both  the  25-  and  1-ps  initial  rise  times,  the  following 
conclusions  seem  plausible.  Not  only  is  the  inertial  spike  a  radial  inertia  effect,  but  so  also  is  the 
reduction  from  this  large  compressive  stress  to  zero  axial  stress,  even  though  the  specimen  is 
subject  to  compressive  axial  strain.  This  reduction  to  zero  axial  stress  is  driven  by  the  rapid  drop 
in  the  axial  acceleration  of  the  incident  bar  to  (and  below)  zero.  Once  the  axial  stress  at  the 
interface  drops  to  zero,42  a  gap  can  form  there.  For  cases  with  more  gradual  deceleration  via 
better  pulse  shaping,  as  in  Moy  et  al.  (5),  there  is  a  less  severe  drop  in  stress  that  never  reaches 
zero  before  the  stress  begins  to  increase  again. 

Finally,  it  should  be  kept  in  mind  that  when  a  gap  exists  over  a  substantial  portion  of  the  face  of 
the  specimen,  the  contact  area  between  the  specimen  and  the  bar  is  reduced,  so  the  local  value  of 
the  axial  stress  in  the  specimen  is  magnified  at  those  points  in  contact  with  the  bar,  whereas  it  is 
zero  at  points  on  the  specimen  face  where  a  gap  has  formed.  Consequently,  in  this  case  the  mean 
axial  stress  computed  as  above  does  not  provide  a  meaningful  estimate  of  the  highly  non-uniform 
stress  state  in  the  specimen.  Furthermore,  differences  between  the  force  at  the  gage  location  and 
the  force  at  the  S-TB  interface  could  arise  due  to  axial  acceleration  of  the  transmission  bar  near 
the  interface  (25).  This  acceleration  is  often  regarded  as  negligible  (based  the  assumption  that 
the  transmitted  pulses  are  weak),  but  it  must  be  significant  at  least  part  of  the  time  in  the  present 


41  The  stress  history  at  the  S-TB  interface  in  Song  et  al.  (6)  was  inferred  from  strain  gages  mounted  on  the  transmission  bar. 
Similarly,  the  stress  history  at  the  S-IB  interface  in  a  SHPB  test  can  be  inferred  from  strain  gages  mounted  on  the  incident  bar. 
Large  stress  oscillations  and  tensile  stresses  are  also  seen  in  some  incident  bar  strain  gage  data  on  soft  materials;  cf.  figure  6 
(lead)  in  Gray  (7)  and  figure  2  (estane-based  polymer  binder)  in  Gray  and  Blumenthal  (5).  This  could  be  indicative  of  gap 
formation  at  the  S-IB  interface.  However,  the  estimate  for  the  stress  at  the  S-IB  interface  based  on  2-wave  calculations  from 
strain  gage  data  is  subject  to  larger  uncertainty  than  the  estimate  for  the  stress  at  the  S-TB  interface;  cf.  Gray  (7),  Gray  and 
Blumenthal  (3),  and  Chen  et  al.  (33).  Quartz  gages  were  used  to  measure  force  histories  in  the  incident  and  transmission  bars  in 
SHPB  tests  on  RTV  630  silicone  rubber  in  Chen  et  al.  (35).  Large  stress  oscillations  and  tensile  stresses  are  seen  in  the  incident 
bar  data  here  as  well  (figures  6  and  8);  and  some  of  the  transmission  bar  data  shows  less  frequent  and  smaller  amplitude  dips  to 
tensile  stress.  These  results  could  be  indicative  of  gap  formation  at  both  interfaces.  However,  the  authors  point  out  that  the 
specimen  may  have  failed  in  some  of  these  tests,  and  this  possibility  complicates  the  interpretation  of  the  data. 

42  Note  that  tensile  axial  stresses  cannot  be  transmitted  across  the  specimen-bar  interface. 
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case.  Indeed,  the  tensile  forces  at  the  gage  location  at  various  times  between  186  and  260  |us 
cannot  possibly  represent  the  force  at  the  S-TB  interface,  since  tensile  forces  cannot  be 
transmitted  across  the  interface:  the  force  on  exerted  on  the  transmission  bar  by  the  specimen 
(and  hence,  on  the  specimen  by  the  transmission  bar)  must  be  compressive  (positive)  or  zero. 


8.  Simulations  with  a  Linear  Elastic  Model  for  the  Specimen 


For  all  of  the  results  presented  up  to  this  point,  the  compressible  Mooney-Rivlin  model  was  used 
for  the  specimen.  Recall  that  this  is  a  nonlinear  elastic  constitutive  model  that  was  calibrated  to 
yield  a  nearly  incompressible  material  (see  section  3  and  appendix  A).  The  question  naturally 
arises  as  to  whether  the  nonlinearity  in  the  material  response  and  the  near  incompressibility  of 
the  material  are  necessary  for  the  formation  of  gaps  or  for  the  persistence  of  the  gaps  once  they 
form.  The  simulations  in  this  section  address  these  issues  by  replacing  the  Mooney-Rivlin  model 
for  the  specimen  with  an  isotropic,  linear  elastic  model  and  by  varying  the  Poisson’s  ratio  for 
that  model.43 

We  used  LS-DYNA’s  Orthotropic  Elastic  model  (Material  Model  2  in  LS-DYNA  [§]  and 
Hallquist  [9])  for  the  specimen,  with  the  material  constants  chosen  in  such  a  way  as  to  yield  an 
isotropic,  linear  elastic  model.  For  all  of  the  simulations  in  this  section,  the  Young’s  modulus  E 
was  fixed  at  0.43  MPa,  which  is  slightly  less  than  twice  the  initial  value  of  the  Young’s  modulus 
used  for  the  Mooney-Rivlin  model  (see  table  1).  The  density  of  the  specimen  remained  the  same 
(1  g/cm  ).  Six  different  values  of  the  Poisson’s  ratio  v  were  used  in  the  simulations,  ranging 
from  0.49999  to  0.48.  The  largest  value,  v  =  0.49999,  is  the  same  value  used  for  the 
compressible  Mooney-Rivlin  model.  The  six  values  for  Poisson’s  ratio  as  well  as  the 
corresponding  values  of  the  other  elastic  constants  determined  from  E  and  v  are  listed  in  table  1. 
Note  that  as  Poisson’s  ratio  decreases  from  0.49999  to  0.48,  the  bulk  modulus  K  and  the 
longitudinal  modulus  L  (which  are  nearly  identical)  decrease  by  a  factor  of  2000,  whereas  the 
changes  in  the  shear  modulus  G  are  insignificant;  the  ratio  of  the  shear  to  bulk  modulus,  G/K , 
which  is  a  measure  of  incompressibility  (see  section  2.3),  increases  from  0.0002  to  0.04. 

8.1  The  1-ps  Initial  Rise  Time 

Figures  31-35  present  results  for  the  loading  wave  with  a  1-ps  initial  rise  time.  Figure  31 
compares  the  histories  of  the  gap  size  on  the  centerline  at  the  two  interfaces  for  the  case 
v  =  0.49999.  Also  plotted  in  that  figure  (right  axis)  is  the  history  of  the  axial  stretch  Az;  see 
equation  8.  Clearly,  nonlinearity  in  the  constitutive  model  is  not  required  for  the  formation  of  a 
gap  or  for  the  persistence  of  the  gap  out  to  large  strains.  The  gap  forms  at  about  165  ps  (3% 
nominal  strain).  The  gap  at  the  S-IB  interface  remains  open  for  a  substantial  time  interval, 


43  The  elastic  constants  used  in  these  simulations  are  not  necessarily  those  appropriate  for  the  small  strain  response  of 
ballistic  gelatin. 
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eventually  closing  at  about  263  ps  (27%  nominal  strain).  The  S-TB  gap  is  also  closed  at  this 
time,  but  has  opened  and  closed  several  times  prior  to  that,  and  the  size  of  the  gap  at  the  S-TB 
interface  during  this  first  stage  is  substantially  less  than  at  the  S-IB  interface.  All  of  these 
features  are  qualitatively  similar  to  those  observed  for  the  nonlinear  model  (see  section  6.1). 


Figure  31.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  and  S-TB  interfaces 
for  a  1-ps  rise  time  and  v  =  0.49999  (linear  elastic  model).  Also  shown  is  the  history  of  mean 
axial  stretch  in  the  specimen  (blue). 
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Figure  32.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for 
a  1-ps  rise  time  and  selected  values  of  Poisson’s  ratio  (linear  elastic  model). 

Figure  32  compares  the  histories  of  the  gap  size  on  the  centerline  at  the  S-IB  interface  for  the  six 
values  of  Poisson’s  ratio.  Figure  33  does  the  same  at  the  S-TB  interface.  Gaps  form  at  both 
interfaces  for  v  >  0.49,  that  is,  for  G/K<  0.02.  But  they  do  not  form  at  either  interface  for  the 
lowest  value  of  Poisson’s  ration,  v  =  0.48,  that  is,  for  G/K  =  0.04.  Thus  it  appears  that  the 
specimen  must  be  nearly  incompressible  in  order  for  gaps  to  form. 
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Figure  33.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for 
a  1-ps  rise  time  and  selected  values  of  Poisson’s  ratio  (linear  elastic  model). 

The  three  values  of  Poisson’s  ratio  that  resulted  in  the  largest  gap  sizes  at  both  interfaces  were 
0.499,  0.4995,  and  0.4999.  Note  that  the  corresponding  gap  histories  in  figures  32  and  33  do  not 
extend  all  the  way  to  360  jus.  We  terminated  these  histories  prior  to  that  time  because  the 
specimen  elements  in  direct  contact  with  the  bar  exhibited  hourglassing.  The  gap  size  curves  in 
the  figures  were  terminated  well  before  this  occurred. 

It  is  interesting  to  note  that  the  maximum  value  (over  the  duration  of  the  simulation)  of  the  size 
of  the  gap  at  the  centerline  does  not  vary  monotonically  with  Poisson’s  ratio.  Figure  34  plots  this 
maximum  gap  size  at  each  interface  as  a  function  of  Poisson’s  ratio  v.  At  the  S-IB  interface,  the 
maximum  gap  size  increases  with  v  up  to  a  maximum  of  140  pm  for  v  =  0.4995,  and  decreases 
with  v  thereafter.  At  the  S-TB  interface,  the  maximum  gap  size  increases  with  v  up  to  a 
maximum  of  180  pm  for  v  =  0.4990,  and  decreases  with  v  thereafter.  For  v  =  0.4995,  the 
maximum  gap  size  at  both  interfaces  is  essentially  the  same  (about  140  pm).  For  smaller  values 
of  v,  the  largest  gap  occurs  at  the  S-TB  interface;  for  larger  values  of  v,  the  largest  gap  occurs  at 
the  S-IB  interface. 
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Figure  34.  The  maximum  gap  size  at  the  centerline  for  each  interface  as  a  function  of  Poisson’s  ratio 
(linear  elastic  model,  1-ps  rise  time). 

The  gap  of  180  pm  on  the  centerline  at  the  S-TB  interface  for  the  case  v  =  0.4990  is  the  largest 
gap  size  observed  in  all  of  the  simulations  in  this  report.  This  occurred  at  276  ps,  which  is  about 
1 10  ps  after  the  loading  wave  arrived  at  the  specimen.  Figure  35  plots  the  deformed  mesh  in  the 
vicinity  of  the  specimen  at  this  instant  for  this  value  of  Poisson’s  ratio.  The  centerline  lies  along 
the  bottom  of  the  plot,  so  the  figure  shows  the  entire  deformed  specimen  mesh  in  the  2-D 
axisymmetric  simulation.  It  is  clear  that  at  this  instant  the  gaps  at  either  interface  are  largest  on 
the  centerline.  Also,  the  gap  extends  further  outward  radially  at  the  S-TB  interface.  Along  most 
of  the  radial  extent  of  the  gap  at  the  S-TB  interface  and  along  half  of  the  radial  extent  of  the  gap 
at  the  S-IB  interface,  the  gap  size  exceeds  the  100- pm  length  of  bar  elements  and  hence  is 
several  times  larger  than  the  50-pm  initial  length  of  the  undeformed  specimen  elements,  which  in 
turn  is  much  larger  than  the  deformed  element  length  at  this  instant.  This  observation  supports 
the  assessment  that  the  gap  phenomenon  cannot  be  attributed  to  numerical  artifact  associated 
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with  an  insufficiently  fine  mesh.  Finally,  note  that  according  to  figure  31,  the  estimated  axial 
stretch  at  this  instant  is  about  0.69;  the  corresponding  nominal  strain  is  31%.  However,  these 
estimates  are  based  on  equations  8  and  9,  that  is,  on  the  distance  between  the  two  bars.  It  is  clear 
from  figure  35  that  for  radial  locations  near  the  centerline,  the  axial  strain  in  the  specimen  is 
substantially  higher  than  this  estimated  value. 


Figure  35.  Deformed  mesh  in  the  vicinity  of  the  specimen  at  t  =  276  ps  for  Poisson’s  ratio 
v  =  0.4990  (linear  elastic  model,  I  -us  rise  time). 
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8.2  The  25- ps  Initial  Rise  Time 


Figures  36-40  present  results  for  the  loading  wave  with  a  25-ps  initial  rise  time.  Figure  36 
compares  the  histories  of  the  gap  size  on  the  centerline  at  the  two  interfaces  for  the  case 
v  =  0.49999.  Also  plotted  in  that  figure  (right  axis)  is  the  history  of  the  axial  stretch  Az.  The  gap 
forms  at  about  185  ps  (5%  nominal  strain).  The  pulse  shaping  has  delayed  the  onset  of  gap 
formation,  has  substantially  reduced  the  gap  size  at  both  interfaces  (the  largest  value  is  now 
about  2.6  pm),  and  has  also  substantially  reduced  the  duration  over  which  the  gap  first  remains 
open  at  the  S-IB  interface.  These  features  are  similar  to  those  observed  for  the  nonlinear  model 
(see  section  7). 
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Figure  36.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  and  S-TB 
interfaces  for  a  25 -ps  rise  time  and  v  =  0.49999  (linear  elastic  model).  Also  shown  is 
the  history  of  mean  axial  stretch  in  the  specimen  (blue). 
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Figure  37.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface 
for  a  25-ps  rise  time  and  selected  values  of  Poisson’s  ratio  (linear  elastic  model). 

Figures  37  compares  the  histories  of  the  gap  size  on  the  centerline  at  the  S-IB  interface  for  the 
six  values  of  Poisson’s  ratio.  Figure  38  does  the  same  at  the  S-TB  interface.  Just  as  for  the  1-ps 
initial  rise  time,  gaps  form  at  both  interfaces  for  v>  0.49,  but  they  do  not  form  at  either  interface 
for  the  lowest  value  of  Poisson’s  ration,  v  =  0.48.  For  simulations  that  exhibited  hourglassing, 
the  gap  history  curves  were  terminated  well  before  hourglassing  occurred. 
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Figure  38.  A  comparison  of  the  histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface 
for  a  25-ps  rise  time  and  selected  values  of  Poisson’s  ratio  (linear  elastic  model). 

Just  as  for  the  1  -jlis  initial  rise  time,  the  maximum  value  (over  the  duration  of  the  simulation)  of 
the  size  of  the  gap  at  the  centerline  does  not  vary  monotonically  with  Poisson’s  ratio.  Figure  39 
plots  this  maximum  gap  size  at  each  interface  as  a  function  of  Poisson’s  ratio  v.  At  both 
interfaces  the  maximum  gap  size  increases  with  v  up  to  v  =  0.499,  and  decreases  with 
v  thereafter.  The  maximum  gap  sizes  for  each  value  of  Poisson’s  ratio  are  roughly  the  same  at 
both  interfaces,  with  the  largest  differences  occurring  for  v  =  0.499. 
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Figure  39.  The  maximum  gap  size  at  the  centerline  for  each  interface  as  a  function  of  Poisson’s 
ratio  (linear  elastic  model,  25-ps  rise  time). 

Figure  40  plots  the  deformed  mesh  in  the  vicinity  of  the  specimen  at  286  ps  for  the  case 
v  =  0.499.  For  this  value  of  v  and  this  instant,  the  gap  size  at  the  S-TB  interface  achieves  its 
largest  value  of  62  pm,  and  the  gap  at  the  S-IB  interface  is  71  pm,  which  is  nearly  equal  to  its 
largest  value  of  72  pm  (see  figures  38  and  39).  From  figure  36,  we  see  that  the  estimated  axial 
stretch  at  this  time  is  about  0.70,  which  implies  a  nominal  strain  of  30%.  The  centerline  lies 
along  the  bottom  of  the  plot,  so  the  figure  shows  the  entire  (deformed)  specimen  mesh  in  the  2-D 
axisymmetric  simulation.  Both  gaps  extend  to  about  half  the  radius  of  the  specimen.  Along 
most  of  the  radial  extent  of  the  gap  at  either  interface,  the  size  of  the  gap  exceeds  the  initial 
length  of  the  specimen  elements  and  is  about  twice  the  deformed  length  of  the  elements. 
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Figure  40.  Deformed  mesh  in  the  vicinity  of  the  specimen  at  t  =  286  ps  for  Poisson’s  ratio 
v  =  0.4990  (linear  elastic  model,  25-ps  rise  time). 
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9.  Discussion  and  Concluding  Remarks 


The  SHPB  test  is  designed  to  impose  a  state  of  compressive  uniaxial  stress  on  the  specimen,  at 
least  after  an  initial  ringing-up  period.  Even  during  the  ring-up,  one  expects  that  the  specimen  is 
in  a  (generally  non-uniform)  state  of  compression,  since  the  initial  compressive  longitudinal 
wave  in  the  specimen  reflects  from  the  higher  impedance  pressure  bars  as  a  compressive  wave. 
Thus  the  implicit  assumption  that  the  specimen  and  pressure  bars  remain  in  contact  would  appear 
to  be  reasonable. 

Nevertheless,  in  this  report  we  have  presented  conclusive  numerical  evidence  that  gaps  may  open 
between  the  specimen  and  the  bars  under  certain  conditions.  These  gaps  formed  at  small  strains 
but  existed  out  to  large  strains,  closing  and  re-opening  multiple  times.  In  some  cases  the  gaps 
persisted  for  over  100  ps  and  extended  over  much  of  a  face  of  the  specimen  for  most  of  that 
time.  The  results  in  the  appendix  support  the  conclusion  that  this  gap  phenomenon  is  not  an 
artifact  of  an  insufficiently  fine  mesh  (appendix  B-l)  or  of  the  contact  algorithm  (appendix  B-2) 
or  of  the  particular  code  used  for  the  simulations  (appendix  C-3).  If  a  gap  does  exist  over  a 
substantial  portion  of  the  specimen  face  in  a  SHPB  test,  then  that  test  is  invalid  for  the  purposes 
of  material  property  characterization  since  the  estimates  for  the  axial  stress,  strain,  and  strain  rate 
are  no  longer  representative  of  the  non-uniform  conditions  in  the  specimen. 

The  results  of  this  study  indicate  that  for  soft  specimens,  gaps  are  most  likely  to  form  under  a 
combination  of  two  conditions,  one  pertaining  to  the  specimen  properties  and  the  other  to  the 
loading  conditions.  For  a  given  loading  condition,  gaps  seem  more  likely  to  form  in  specimens 
which  are  nearly  incompressible.  This  condition  was  satisfied  in  the  two  simulations  with  the 
nonlinear  elastic  (Mooney-Rivlin)  model  for  the  specimen,  and  gaps  formed  in  both  cases 
(sections  6  and  7).  When  a  linear  elastic  model  was  used  for  the  specimen,  gaps  formed  for  five 
of  the  six  sets  of  elastic  moduli  considered  (section  8).  The  case  where  no  gaps  formed 
corresponded  to  the  lowest  Poisson’s  ratio  (0.48)  and  hence  to  the  highest  ratio  of  shear  to  bulk 
modulus  (0.04). 

The  loading  condition  that  seems  to  promote  gap  formation  is  axial  deceleration  of  the  specimen. 
As  the  time  rate  of  change  of  the  velocity  imposed  by  the  incident  bar  on  the  specimen  decreases 
from  peak  positive  values  (acceleration)  through  zero  to  negatives  values  (deceleration),  the  axial 
stress  at  the  specimen-bar  interfaces  drops  to  zero,  allowing  the  specimen  to  separate  from  the 
bars.  It  appears  that  the  driving  mechanism  behind  this  inertial  effect  is  the  corresponding  radial 
deceleration  induced  by  the  axial  deceleration.  This  conclusion  is  not  inconsistent  with  existing 
literature  on  radial  inertial  effects  for  soft  specimens  (26-33).  Furthermore,  we  found  that 
constraining  the  lateral  surface  of  the  specimen  (to  eliminate  the  radial  motion)  suppresses  gap 
formation,  at  least  in  direct  impact  tests  (appendix  C-l.l). 
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On  the  other  hand,  the  simple  qualitative  conclusions  above  belie  the  extremely  complex  stress 
states  in  the  specimen  and  pressure  bars  (particularly  the  incident  bar)  at  times  in  the  vicinity  of 
gap  formation,  as  illustrated  by  the  pressure  contour  plots  in  section  6.3.44  Starting  about  1.6  ps 
prior  to  the  opening  of  gaps  at  both  interfaces  and  at  least  1  ps  prior  to  the  onset  of  negative 
pressure  in  the  specimen,  the  pressure  in  the  incident  bar  oscillates  (with  axial  distance  from  the 
specimen)  between  positive  and  negative  values,  eventually  becoming  negative  in  the  entire 
vicinity  of  the  specimen  immediately  before  the  gaps  form  (section  6.3). 

To  the  best  of  our  knowledge,  the  gap  phenomenon  documented  here  has  not  been  reported  in 
either  the  experimental  or  the  computational  literature.  However,  certain  features  of  previously 
reported  experimental  data  appear  to  be  consistent  with  the  opening  and  closing  of  gaps,  namely, 
axial  stress  histories  (inferred  from  strain  or  stress  gage  measurements  in  the  bars)  which 
oscillate  between  positive  and  negative  values  (see  section  7.3).  Since  our  initial  concern  was  to 
determine  whether  or  not  these  gaps  were  a  numerical  artifact,  we  focused  on  examining  a  few 
cases  in  detail  rather  than  simulating  a  broad  range  of  loading  conditions  and  constitutive  models 
for  the  specimen.  Consequently,  our  study  is  by  no  means  exhaustive.  We  conclude  by  listing 
some  of  the  limitations  of  this  work,  which  in  turn  suggest  areas  for  further  study: 

1)  We  certainly  expect  that  gaps  will  no  longer  form  if  the  strain  rate  in  the  specimen  is 
sufficiently  small,  but  only  one  axial  strain  rate  was  considered  here,  a  nominal  rate  of 
2500/s  (after  the  initial  ring-up).  The  imposed  strain  rate  can  be  reduced  either  by 
decreasing  the  plateau  velocity  imposed  at  the  far  end  of  the  incident  bar  or  by  increasing 
the  thickness  of  the  specimen. 

2)  We  considered  only  a  linear  ramp  to  a  plateau  velocity  in  the  main  body  of  the  report  (and 
more  severe  loadings  in  the  appendices  C-l  and  C-3).  It  would  be  of  interest  to  consider 
smoother  loading  histories  at  the  far  end  of  the  incident  bar,  since  experimental  evidence 
indicates  that  gaps  do  not  form  with  better  pulse  shaping. 

3)  We  considered  only  one  specimen  diameter  (12.7  mm).  Theories  on  radial  inertial  effects 
in  SHPB  tests  predict  that  the  inertial  stresses  are  proportional  to  the  square  of  the 
specimen  diameter.  Thus  if  radial  inertia  is  a  mechanism  for  gap  formation,  we  would 
expect  that  gaps  would  no  longer  form  if  the  specimen  diameter  is  sufficiently  small. 

4)  All  the  simulations  in  this  report  involved  specimens  that  were  “soft”  relative  to  the  bars 
(particularly  in  shear),  so  we  cannot  draw  any  conclusions  about  gap  formation  in  stiffer 
specimens. 


44  Recall  that  these  pressure  contours  were  for  the  loading  wave  with  a  1-ps  initial  rise  time. 
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5)  Only  elastic  (linear  or  nonlinear)  constitutive  models  were  considered  for  the  specimen, 
since  absence  of  strain  rate  effects  in  the  model  allowed  for  easier  interpretation  of  the 
results.  It  is  not  clear  how  the  incorporation  of  viscoelasticity  in  the  model  would  affect 
the  formation  of  gaps,  although  we  have  no  reason  to  believe  that  this  would  suppress  gap 
formation  altogether. 

6)  For  the  nonlinear  elastic  specimen  model  (the  compressible  Mooney-Rivlin  model),  we 
considered  only  one  set  of  material  parameters.  It  would  be  instructive  to  fix  the  bulk 
modulus  and  increase  the  shear  modulus  until  gap  formation  is  suppressed.  Note  that  in  the 
simulations  with  a  linear  elastic  specimen,  we  essentially  fixed  the  shear  modulus  and 
decreased  the  bulk  modulus  until  gap  formation  was  suppressed. 

7)  As  an  approximation  to  a  well  lubricated  specimen,  the  specimen-bar  interfaces  were 
treated  as  frictionless  in  the  simulations.  It  is  not  clear  what  effect  (if  any)  a  small  amount 
of  friction  would  have  on  the  formation  of  gaps. 
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Appendix  A.  The  Mooney-Rivlin  Constitutive  Model 


A-l  Theoretical  Background 

This  subsection  contains  a  brief  discussion  of  the  concepts  and  notation  needed  to  describe  the 
two  nonlinear  elastic  models  considered  below;  for  additional  background  the  reader  may  refer 
to  the  books  (16-19). 

Let  A],  A2,  A3  denote  the  principal  stretches,  that  is,  the  ratios  of  the  deformed  to  undeformed 
length  along  the  principal  axes  of  strain.  Each  stretch  A,,  is  unity  in  the  undeformed  state,  with 
Aj  >  1  in  tension,  and  0  <  Aj  <  1  in  compression.  The  principal  nominal  (engineering)  strains  are 
A— 1  when  taken  positive  in  tension,  or  1  -A,-  when  taken  positive  in  compression.  The  left 
Cauchy-Green  deformation  tensor  B  is  defined  in  terms  of  the  deformation  gradient  F  by 
B  =  FF  ,  where  the  superscript  T  denotes  the  transpose.  The  principal  axes  of  B  are  the 

principal  axes  of  strain  in  the  deformed  state.  The  principal  values  (or  eigenvalues)  of  B  are 

2  2  2 

A]  ,  A2  ,  A3  ,  the  squares  of  the  principal  stretches.  The  three  principal  invariants  of  B  are 

Ix  =  txB  =  A^  +A^+A^,  (A-l) 

I2  =  | [( tr  Bf  -  tr  B 2  ]  =  A{A^  +  A; A;  +  A; A,2 ,  (A-2) 

I3  =  det  B=  AfAfA,2  =  J2.  (A- 3) 

Here  “tr”  denotes  the  trace,  “det”  denotes  the  determinant,  and 

J  =  det  F  =  A 1 A2A3  (A-4) 

is  the  Jacobian  of  the  deformation,  that  is,  the  local  ratio  of  deformed  to  undeformed  volume. 
Thus,  /  =  1  at  points  where  there  is  no  volume  change.  A  useful  alternative  expression  for  L  is 

/,_  =  /,  1  1|  (A-5) 


The  principal  stresses  (i.e.,  the  eigenvalues)  of  the  Cauchy  stress  tensor  cr  are  denoted  by 
ai,  02,  cr3,  and  the  pressure  p  is  given  by 


P  = 


1  1/ 

3trCT  =  ~r1+G2 


+  g3). 


(A-6) 


For  any  (not  necessarily  linear)  isotropic  elastic  material,  the  principal  axes  of  stress  coincide 
with  the  principal  axes  of  B. 
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A-2  The  Mooney-Rivlin  Model  for  Incompressible  Elastic  Materials 

The  classical  Mooney-Rivlin  model  is  a  nonlinear  elastic  constitutive  model  for  isotropic, 
incompressible  solids;  cf.  the  original  papers  by  Mooney  (13)  and  Rivlin  ( 14),  as  well  as  the 
books  (15-19).  Since  the  material  is  assumed  to  be  incompressible,  there  is  no  change  in 
volume,  so  equations  A-3  and  A-4  reduce  to 

I3  =  J  =  X]Xih=\.  (A-7) 

The  strain  energy  function  W  for  the  Mooney-Rivlin  model  has  the  simple  form 

W  =  A1(/1-3)  +  A,(/2-3),  (A- 8) 

The  coefficients  A  i  >  0  and  At  >  0  are  constants  with  the  units  of  stress  (i.e.,  elastic  moduli);  the 
inequalities  are  required  for  physically  realistic  behavior.  The  Cauchy  stress  tensor 
corresponding  to  this  strain  energy  function  is 

ct  =  —p  I  +  2AlB  -  2A2B~l ,  (A-9) 


where  I  denotes  the  identity  tensor.  The  shear  modulus  G  in  the  small  strain,  linear  elastic 
approximation  to  the  constitutive  relation  A-9  is  given  by  G  =  2(A  i  +  At).  The  special  case  of 
equation  A-9  with  At  =  0  (and  hence  2Ai  =  G  )  is  called  a  neo-Hookean  material. 

The  coefficient  p  in  equation  A-9  denotes  an  indeterminate  scalar  related  to  the  pressure.  The 
indeterminacy  of  p  is  a  consequence  of  the  incompressibility  constraint:  p  is  not  determined 
by  B ,  but  the  value  of  p  at  each  place  and  time  must  be  such  that  the  momentum  balance 
equations  and  the  boundary  conditions  are  satisfied.  Although  p  in  equation  A-9  is  often 
referred  to  as  a  “pressure”,  it  is  generally  not  equal  to  the  pressure  p  as  defined  in  equation  A-6, 
since  the  other  terms  in  equation  A-9  contribute  to  the  trace  of  a. 


The  relation  A-9  implies  that  the  principal  stresses  are  given  by 

l 


oi=-p  +  2AlAf-2A2  —  (i  =  1,2,3). 


The  difference  of  any  two  principal  stresses  is  independent  of  p  : 

o;-cyj=2A1(A>-Aj)-2A2 


'l  r 


(A- 10) 


(A- 11) 


Uniaxial  Stress:  The  general  relation  A-l  1  may  be  used  to  obtain  the  stress-stretch  relation  in  a 
uniaxial  stress  test.  Consider  a  cylindrical  coordinate  system  with  the  direction  of  applied  stress 
parallel  to  the  z-axis.  Then  azz  is  the  only  non-zero  stress  component.  If  azz  >  0,  this  is  a  simple 
tension  test;  if  ozz  <  0,  it  is  a  simple  compression  test.  Since  the  material  is  isotropic,  the 
principal  axes  of  strain  lie  along  the  coordinate  axes  in  this  case  and  the  strains  (and  hence 
stretches)  in  the  radial  and  hoop  directions  are  equal.  Thus,  in  the  relations  above,  we  may  set 
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A\  =  Ar  =  Ae  =  A2  and  A3  =  Xz,  as  well  as  ai  =  o,T  =  aee  =  a2  =  0  and  as  =  aK.  Then  by  equations 
A-6  and  A-7,  we  see  that  p  =  -azz/ 3  and 


A„  —  Aa  — 


41: 


(A- 12) 


Then  on  setting  i  =  3  and  j  =  1  in  equation  A-l  1  and  using  the  relations  above,  we  obtain  the 
following  expression  for  the  axial  stress  azz  in  terms  of  the  axial  stretch  Az: 


f 


=  2A 


1 


At - 

v  Z  Ay 


-24 


3G(A- 1). 


(A- 13) 


The  expression  on  the  right  is  the  linearly  elastic  approximation,  which  is  valid  for  small  strains 
only,  that  is,  Az  close  to  1.  For  an  incompressible  material,  3 G  =  E,  where  E  is  the  Young’s 
modulus,  so  we  recover  the  well-known  result  that  in  the  linear  elastic  approximation,  the  axial 
stress  is  the  Young’s  modulus  times  the  axial  strain. 


Stress  Difference  on  the  Centerline:  Consider  an  axisymmetric  (but  possibly  non-uniform) 
deformation  of  a  solid  cylindrical  specimen.  As  discussed  at  the  end  of  section  2.6,  at  points  on 
the  centerline  (i.e.,  the  z-axis),  the  radial,  hoop  and  axial  directions  are  principal  axes  of  stress 
and  strain,  and  the  relations  14  in  section  2.6  hold  exactly.  Thus,  the  principal  stress  difference 
ozz  -  <jrr  is  given  by  equation  A-l  1  with  the  indices  i  and  j  set  to  z  and  r,  respectively;  and  since 
7=1,  equation  A- 12  holds.  Therefore, 


f 


azz  _  °rr  ~  2A, 


1 


At - 

v  Z  Ay 


-2A, 


3G(A  -1) 


(A- 14) 


where  the  linearly  elastic  approximation  on  the  right  holds  for  small  strains. 

Note  that  the  right-hand  side  of  this  equation  is  the  same  as  that  of  equation  A- 13.  In  particular, 
for  a  uniaxial  stress  test,  equation  A-14  reduces  to  equation  A-13  since  <jrr  =  0  in  that  case. 
However,  for  the  more  general  case  considered  here,  o„-  need  not  be  zero,  and  equation  A-14  is 
only  valid  on  the  centerline. 

A-3  The  Compressible  Version  of  the  Mooney-Rivlin  Model  in  LS-DYNA 

The  simulations  described  in  sections  6  and  7  used  LS-DYNA’s  compressible  version  of  the 
Mooney-Rivlin  model  for  the  specimen;  this  is  Material  Type  27  in  (5,  9).  The  model  has  three 
material  constants  (aside  from  the  density):  the  elastic  moduli45  A\  and  and  the  Poisson’s  ratio 
v.  The  strain  energy  (per  unit  undeformed  volume)  is  given  by 

W  =  A1(/1-3)  +  A2(/2-3)  +  C(/3' 2  — l)  +  Z)(/3  — l)2 ,  (A- 15) 


45  The  LS-DYNA  manuals  (8,  9)  use  the  symbols  A  and  B  for  the  constants  A{  and  A2,  respectively. 
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which  differs  from  the  incompressible  model  by  the  addition  of  the  terms  with  coefficients  C  and 
D.  The  constant  C  is  determined  from  A  i  and  A2  by 

C  =  ^-  +  A2.  (A- 16) 

The  constant  D  is  determined  from  A\,  A2,  and  v  by 

A,(5v-2)  +  A,(llv-5) 

D  =  — - - - - — - -.  (A- 17) 

2(1 -2v) 

We  refer  to  this  as  the  “compressible  Mooney-Rivlin  model”. 

The  Cauchy  stress  tensor  corresponding  to  this  strain  energy  function  is46 

a  =  -pl  +  2^-B-2A1JB-' ,  (A- 18) 

where 

p=4D(y-/J)  +  i£-^/,.  (A-19) 

When  7=1,  equation  1 8  yields  the  relations  A.  10  and  A.ll  for  the  principal  stresses,  which  were 
obtained  previously  for  the  incompressible  model.  If  7  is  close  to  1,  then  these  relations  hold  as 
approximations.  Similarly,  when  7  =  1,  we  recover  the  relation  A- 13  for  the  axial  stress  azz  in  a 
uniaxial  stress  test  and  the  relation  A- 14  for  the  stress  difference  ozz  -  a,T  on  the  centerline  for  a 
general  axisymmetric  deformation.  If  7  is  close  to  1,  then47  Ar  « 1/  JX" ,  and  we  find  that  A- 13 

and  A- 14  hold  as  approximations.  In  this  regard,  note  that  3 G  ~  E  for  a  nearly  incompressible 
material;  indeed,  for  the  ballistic  gelatin  specimen  they  agree  to  at  least  three  significant  figures 
(see  table  1). 

Note  that  the  scalar  p  in  equations  A-18-A-19  is  generally  not  equal  to  the  pressure  p,  since  the 
other  terms  in  equation  A- 18  contribute  to  the  trace  of  a.  Indeed,  from  equations  A-6,  A- 18, 

A-5,  and  A-19,  we  obtain 

p  =  p--  —  Il+~^I2.  (A-20) 

3  7  3  7 


4^  None  of  the  relations  in  this  paragraph  appear  in  the  LS-DYNA  User’s  Manual  (8)  or  in  the  Theory  Manual  by  Hallquist 
(9).  However,  equation  A- 18  is  equivalent  to  the  relation  for  the  second  Piola-Kirchhoff  stress  tensor  given  in  (9),  and  equations 
A-19-A-20  imply  the  relation  for  the  pressure  in  a  state  of  pure  dilation  given  in  (9). 

47  See  equation  15  in  section  2.6. 
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Thus  the  pressure  p  is  a  complicated  function  of  both  the  volumetric  strain  and  the  shear  strain. 

In  particular,  p  does  not  necessarily  reduce  to  zero  when  there  is  no  volume  change  (i.e.,  when 
7  =  1).  However,  p  does  equal  zero  in  the  undeformed  state,  as  expected.  This  follows  from 
equations  A-20,  A- 19,  A- 16,  and  the  fact  that  7=1  and  I\  =  h  =  3  when  the  material  is 
undeformed. 

For  nearly  incompressible  materials  (v  close  to  Vi),  A\  and  A2  (and  hence  C)  are  small  compared 
to  D,  so  for  small  volume  changes  (7  close  to  but  not  equal  to  1),  equations  A- 19  and  A-20  yield 
the  simple  approximations 

p  «  p  » K(l- 7) ,  K  « 8D ,  (A-21) 

where  the  constant  K  is  the  initial  (or  linear  elastic)  bulk  modulus.  It  is  easily  verified  that  the 
approximation  for  K  in  equation  A-21  is  consistent  with  equations  2i,  16,  and  A- 17  when  v  is 
close  to  Vi. 

The  LS-DYNA  manuals  do  not  provide  any  explanation  or  references  for  the  unconventional 
form  of  the  two  compressibility  terms  in  the  strain  energy  function  in  equation  A- 15,  that  is,  the 
terms  with  coefficients  C  and  D ,  although  the  term  involving  C  appears  to  have  been  added  so 
that  the  pressure  reduces  to  zero  in  the  undeformed  state.  However,  the  relation  for  the  pressure 
obtained  from  this  strain  energy  function  (equations  A-19-A-20)  is  inappropriate  for  large 
volumetric  strains,  since  the  7  -  7 3  term  in  equation  A- 19  has  an  absolute  maximum  at 

7  =  -J=  «  0.58  .  (A-22) 

s 

If  v  is  sufficiently  close  to  Vi  (e.g.,  v  >  0.4999),  the  7  -73  term  dominates  except  at  very  large 
compressions  (i.e.,  7  <  0.5).48  This  physically  unrealistic  behavior  was  not  an  issue  for  the 
SHPB  simulations  reported  here.  For  the  parameters  used  in  our  calibration  of  the  model  (see 
table  3),  the  value  of  the  pressure  at  the  maximum  described  above  is  about  770  MPa,  whereas 
the  peak  pressures  in  our  simulations  were  around  4  MPa,  with  the  exception  of  the  sharp  spikes 
on  the  centerline  which  did  not  exceed  10  MPa.49  In  other  words,  the  volume  changes  were 
small  enough  that  the  approximation  A-21  was  valid. 

A-4  Verification  of  the  Compressible  Mooney-Rivlin  Model 

Prior  to  our  simulations  of  SHPB  tests,  we  performed  a  partial  verification  of  the  implementation 
of  the  compressible  version  of  the  Mooney-Rivlin  model  in  LS-DYNA.  We  used  the  fact  that 
for  a  quasi-static  uniaxial  stress  test,  the  stress  predicted  by  the  compressible  model  should 
approach  that  predicted  by  the  incompressible  model  as  Poisson’s  ratio  v  approaches  Vi  with 
A  i  and  A 2  fixed.  The  constants  Aj  and  A 2  were  chosen  as  in  table  3.  We  used  a  single  8-node 


48  This  unrealistic  behavior  is  caused  by  the  last  term  in  equation  A. 15;  it  could  be  eliminated  by  replacing  I3  there  by  its 
square  root,  that  is.  by  J,  but  then  the  term  with  coefficient  C  would  also  need  to  be  modified  to  guarantee  that  the  pressure  is 
zero  in  the  undeformed  state. 

49  Refer  to  figures  11.  12,  17-23,  27,  and  28,  and  recall  that  in  the  axial  stress  history  plots,  the  stress  state  at  the  stress  peaks 
prior  to  the  formation  of  gaps  is  nearly  hydrostatic. 
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hexagon  element  that  was  initially  cubic  with  a  100  |um  edge  length.  In  order  to  diminish  any 
inertial  effects,  the  initial  density  was  reduced  from  1  to  0.001  g/cm  for  these  simulations  only. 
There  were  no  discemable  differences  in  the  results  for  nominal  axial  strain  rates  of  100/s, 

1000/s,  and  10,000/s. 

For  various  values  of  Poisson’s  ratio  v  between  0.495  and  0.49999,  the  stress  components  were 
computed  for  axial  stretches  Az  ranging  from  0.2  (large  compression)  to  4.0  (large  extension). 
These  results  were  compared  with  the  theoretical  relation  A- 13  for  the  incompressible  model. 

The  results  for  v  =  0.495  are  shown  in  figure  A-l.  The  stress-stretch  curves  for  the  axial  stress 
azz  for  the  compressible  model  (red)  and  incompressible  model  (dashed  line)  are 
indistinguishable  except  for  a  slight  difference  when  Az  is  below  0.35.  This  difference  decreased 
for  higher  values  of  Poisson’s  ratio;  in  particular,  for  the  value  v  =  0.49999  used  in  the  SHPB 
simulations,  the  stress -stretch  curves  were  indistinguishable  for  all  values  of  Az.  As  expected,  the 
other  stress  components50  (the  blue  line  in  the  figure)  were  indistinguishable  from  zero  on  this 
scale. 


Figure  A-l.  Stress  components  as  a  function  of  axial  stretch  A.  for  a  uniaxial  stress 
test  on  a  single  element  using  the  compressible  Mooney-Rivlin  model 
with  v  =  0.495  (colored  curves).  Comparison  with  axial  stress  azz  for  the 
incompressible  model  (dashed  curve). 

Finally,  note  that  for  small  volume  changes  ( /  « 1 ),  the  curve  in  figure  A-l  also  represents  the 
stress  difference  ozz-  <jrr  on  the  centerline  for  a  general  axisymmetric  deformation.51 


50  These  components  are  taken  relative  to  a  Cartesian  coordinate  system  aligned  with  the  (initially  cubic)  specimen. 

51  See  equation  A- 14  and  the  discussion  following  equation  A- 19. 
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Appendix  B.  Sensitivity  Studies 


All  of  the  results  presented  in  the  main  body  of  this  report  were  obtained  with  the  baseline 
25x50-pm  mesh  in  the  specimen;  see  section  4.1  and  figure  2a.  All  of  the  results  presented  in 
the  main  body  were  obtained  with  the  LS-DYNA  default  values  of  the  contact  algorithm 
parameters  (SFAC  =  1,  VDC  =  10);  see  section  4.2.  Here  we  examine  the  sensitivity  of  the  gap 
size  on  the  centerline  to  the  size  of  the  specimen  mesh  (appendix  B-l)  and  the  contact  algorithm 
parameters  (appendix  B-2). 

For  the  mesh  sensitivity  study  the  default  values  of  the  contact  algorithm  parameters  are  used. 
Likewise,  for  the  contact  algorithm  sensitivity  study,  the  baseline  mesh  is  used.  For  all  of  the 
results  in  this  appendix  we  used  the  compressible  Mooney-Rivlin  constitutive  model  for  the 
specimen  and  the  larger  (12. 8 -mm)  bar  radius.  Both  the  1-  and  25 -ps  initial  rise  times  for  the 
loading  wave  are  considered. 

B-l  Mesh  Sensitivity 

Refer  to  section  4.1  for  a  discussion  of  the  specimen  and  bar  meshes.  Figures  B-l-B-4  compare 
the  histories  of  the  gap  size  on  the  centerline  for  three  specimen  meshes:  the  baseline 
25x50  pm  mesh  (figure  2a),  the  coarser  50x100pm  mesh  (figure  2b),  and  the  finer 
12.3x12.5  pm  mesh  (figure  2c).  The  mesh  size  in  the  bars  is  100x100  pm  for  all  simulations. 
Figures  B-l  and  B-2  are  for  the  1-ps  initial  rise  time  for  the  loading  wave;  figures  B-3  and  B-4 
are  for  the  25-ps  initial  rise  time.  Figures  B-l  and  B-3  give  the  gap  size  at  the  S-IB  interface; 
figures  B-2  and  B-4  give  the  gap  size  at  the  S-TB  interface. 


85 


Figure  B-l. 


Histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for  three  specimen  meshes 
and  a  1-ps  initial  rise  time  for  the  loading  wave. 
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Figure  B-2.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for  three  specimen 
meshes  and  a  1  -us  initial  rise  time  for  the  loading  wave. 
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Figure  B-3. 


Histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for  three  specimen  meshes 
and  a  25-ps  initial  rise  time  for  the  loading  wave. 
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Figure  B-4.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for  three  specimen  meshes 
and  a  25-ps  initial  rise  time  for  the  loading  wave. 

First,  we  consider  the  1-jlis  initial  rise  time  for  the  loading  wave.  A  detailed  discussion  of  the 
results  for  this  case  using  the  baseline  mesh  was  given  in  section  6.1,  where  the  gap  formation 
was  broken  down  into  six  stages.  At  both  interfaces,  there  is  mesh  convergence  in  the  onset  of 
the  gaps  (at  165  ps).  For  the  gap  at  the  S-IB  interface  (figure  B-l),  we  see  that  there  is 
substantial  mesh  convergence  in  both  gap  size  and  duration  for  Stage  II  (165-280  ps,  during 
which  there  is  a  persistent  gap)  and  Stage  III  (280-286  ps,  during  which  the  gap  remains  closed). 
For  the  gap  at  the  S-TB  interface  (figure  B-2),  we  see  that  there  is  substantial  mesh  convergence 
in  the  duration  of  Stages  II  and  III,  but  no  convergence  in  the  size  of  the  gap  during  Stage  II. 
However,  the  results  for  all  three  meshes  are  qualitatively  similar  during  Stage  II,  that  is,  the  gap 
opens  and  closes  several  times.  For  Stage  IV  (286-324  ps,  during  which  both  gaps  re-open), 
there  is  a  degree  of  convergence  in  the  duration  of  the  gaps  at  either  interface  but  not  in  gap  size. 
Convergence  has  not  been  achieved  at  later  times. 
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Next,  we  consider  the  25-ps  initial  rise  time  for  the  loading  wave.  A  discussion  of  the  results  for 
this  case  using  the  baseline  mesh  was  given  in  section  7.  At  both  interfaces,  there  is  mesh 
convergence  in  the  onset  of  the  gaps.  For  the  gap  at  the  S-IB  interface  (figure  B-3),  we  see  that 
there  is  substantial  mesh  convergence  in  both  gap  size  and  duration  during  the  Stage  II 
(184-218  ps,  during  which  there  is  a  persistent  gap);  in  particular,  the  timing  of  the  peak  value 
(2. 2-2.4  pm)  agrees  within  1  ps  for  all  three  meshes.  For  other  times  at  the  S-IB  interface  and 
all  times  at  the  S-TB  interface  (figure  B-4),  the  gap  sizes  remain  under  about  1.2  pm;  the  local 
peaks  in  gap  size  generally  occur  at  similar  times;  and  the  gaps  persist  for  similar  durations  for 
the  three  meshes.  However,  gap  sizes  generally  display  less  mesh  convergence  than  gap 
durations. 

The  results  above  support  the  conclusion  that  the  gap  phenomenon  is  not  a  numerical  artifact 
introduced  by  an  insufficiently  fine  mesh  in  the  specimen. 

B-2  Sensitivity  to  the  Contact  Algorithm  Parameters 

Refer  to  section  4.2  for  a  brief  discussion  of  the  contact  algorithm.  Figures  B-5,  B-6,  B-9,  and 
B-10  compare  the  histories  of  the  gap  size  on  the  centerline  for  three  values  of  the  scale  factor 
SFACT:  0.1,  1  and  10,  which  correspond  to  a  less  stiff,  the  default,  and  a  stiffer  penalty  force, 
respectively.  Figures  B-7,  B-8,  B-ll,  and  B-12  compare  the  histories  of  the  gap  size  on  the 
centerline  for  two  values  of  the  viscous  damping  coefficient  VDC:  10  and  20,  which  correspond 
to  the  default  damping  and  twice  the  default  damping,  respectively.  Figures  B-5-B-8  are  for  the 
1-ps  initial  rise  time  for  the  loading  wave;  figures  B-9-B-12  are  for  the  25-ps  initial  rise  time. 
The  results  at  the  S-IB  interface  are  given  in  figures  B-5,  B-7,  B-9,  and  B-ll;  those  for  the  S-TB 
interface  are  given  in  figures  B-6,  B-8,  B-10,  and  B-12. 

First,  we  consider  the  1-ps  initial  rise  time  for  the  loading  wave  (figures  B-5-B-8).  A  detailed 
discussion  of  the  results  for  this  case  using  the  default  values  for  the  contact  parameters  was 
given  in  section  6.1,  where  the  gap  formation  was  broken  down  into  six  stages.  The  durations  of 
Stage  II  (165-280  ps),  Stage  III  (280-286  ps)  and  Stage  IV  (286-324  ps)  are  unaffected  by  the 
changes  in  SFACT  and  VDC.  The  size  of  the  gap  at  the  S-IB  interface  during  Stage  II 
(165-280  ps)  is  also  unaffected  by  the  changes  in  SFACT  and  VDC.  At  other  times,  the  size  of 
the  gap  at  both  interfaces  varies  somewhat  with  the  values  of  these  parameters,  but  the  results  are 
qualitatively  similar. 


90 


Figure  B-5.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for  three  values  of  the 
contact  parameter  SFACT  and  a  1-ps  initial  rise  time  for  the  loading  wave. 


91 


Figure  B-6.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for  three  values  of  the 
contact  parameter  SFACT  and  a  1-ps  initial  rise  time  for  the  loading  wave. 
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Figure  B-7.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for  two  values  of  the 
contact  parameter  VDC  and  a  I  -us  initial  rise  time  for  the  loading  wave. 
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Figure  B-8.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for  two  values  of  the 
contact  parameter  VDC  and  a  I  -us  initial  rise  time  for  the  loading  wave. 

Next,  we  consider  the  25-ps  initial  rise  time  for  the  loading  wave  (figures  B-9-B-12).  A 
discussion  of  the  results  for  this  case  using  the  default  values  for  the  contact  parameters  was 
given  in  section  7.  At  the  S-IB  interface,  we  see  that  doubling  VDC  has  no  effect  on  the  duration 
or  size  of  the  gap  during  the  Stage  II  (184-218  ps).  The  changes  in  SFACT  also  have  no  effect 
on  the  duration  and  only  a  slight  effect  on  the  size  of  the  gap  during  this  stage;  in  particular,  the 
peak  value  occurs  at  about  200  ps  for  all  three  cases.  For  other  times  at  the  S-IB  interface  and  all 
times  at  the  S-TB  interface,  doubling  VDC  has  little  effect  on  the  duration  of  the  gaps  and  only  a 
slight  effect  on  their  size.  The  changes  in  SFACT  have  more  effect  on  both  duration  and  size, 
but  the  results  are  qualitatively  similar. 

The  results  above  support  the  conclusion  that  the  gap  phenomenon  is  not  a  numerical  artifact 
introduced  by  the  contact  algorithm. 
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Figure  B-9.  Histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for  three  values  of  the 
contact  parameter  SFACT  and  a  25-us  initial  rise  time  for  the  loading  wave. 
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Figure  B-10. 


Histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for  three  values  of  the 
contact  parameter  SFACT  and  a  25-ps  initial  rise  time  for  the  loading  wave. 
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Figure  B-ll. 


Histories  of  the  gap  sizes  on  the  centerline  at  the  S-IB  interface  for  two  values  of  the 
contact  parameter  VDC  and  a  25 -us  initial  rise  time  for  the  loading  wave. 
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Figure  B-12. 


Histories  of  the  gap  sizes  on  the  centerline  at  the  S-TB  interface  for  two  values  of  the 
contact  parameter  VDC  and  a  25 -us  initial  rise  time  for  the  loading  wave. 
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Appendix  C.  Additional  Computational  Studies 


This  section  contains  summaries  of  additional  numerical  simulations  in  which  gaps  formed  at  the 
specimen-bar  interfaces.  These  simulations  involve  conditions  that  differ  in  one  or  more  ways 
from  those  considered  previously,  and  in  appendix  C-3  they  involve  the  use  of  a  different  code. 

C-l  Direct  Impact  on  the  Specimen 

Recall  that  in  a  conventional  SHPB  test,  the  striker  bar  impacts  the  incident  bar  and  generates  a 
compressive  stress  pulse  which  propagates  along  the  incident  bar  and  subsequently  into  the 
specimen  and  the  transmission  bar;  and  for  soft  specimens,  a  pulse  shaper  is  often  inserted 
between  the  striker  and  incident  bars  to  produce  a  smoother  loading  wave.  A  variation  of  the 
SHPB  test,  known  as  the  direct  impact  Hopkinson  bar  or  direct  impact  compression  test,  involves 
the  elimination  of  the  incident  bar,  so  that  the  striker  bar  directly  impacts  the  specimen  (27-30). 
The  transmission  bar  is  retained  and  serves  the  same  purpose  as  in  a  conventional  SHPB  test. 

The  main  advantage  of  this  technique  is  that  it  permits  higher  strain  rates  to  be  imposed  on  the 
specimen.  There  are  some  disadvantages  as  well;  see  the  papers  cited  above  and  also  the 
discussions  in  Gray  ( 1 )  and  Jia  and  Ramesh  (31). 

We  performed  several  simulations  in  which  a  striker  bar  directly  impacts  the  specimen.  The 
specimen  geometry,  bar  diameters,  bar  properties,  and  mesh  sizes  were  the  default  values  used 
previously  (sections  2.1,  2.3,  and  4.1).  The  compressible  Mooney-Rivlin  model  (section  3)  was 
used  for  the  specimen.  The  initial  velocity  of  the  striker  bar  is  vsb  =  3.625  m/s,  which  yields  a 
nominal  axial  strain  rate  of  approximately  vsb/Ts  =  2500/s,  as  before.  Note  that  vSb  is  twice  the 
plateau  velocity  Vo  imposed  at  the  far  end  of  the  incident  bar  in  the  previous  SHPB  simulations; 
on  the  other  hand,  the  particle  velocity  in  the  incident  bar  nearly  doubles  as  the  loading  wave 
reflects  from  the  S-IB  interface  (section  5.3).  Thus,  in  either  case  an  axial  velocity  of  roughly 
3.625  m/s  is  ultimately  imposed  on  the  face  of  the  specimen.  However,  in  the  SHPB  simulations 
the  imposed  velocity  ramps  up  smoothly  from  zero  (figures  5  and  6),  whereas  the  direct  impact 
imposes  an  initial  step  in  velocity  (and  hence  an  acceleration  impulse)  on  the  specimen.  In  the 
subsequent  discussion,  the  SHPB  simulation  mentioned  for  purposes  of  comparison  will  be  the 
case  with  a  1-ps  initial  rise  time  (section  6),  as  this  case  imposed  the  most  severe  accelerations 
on  the  specimen. 

In  the  SHPB  simulation,  the  velocity  at  the  S-IB  interface  exhibits  large  oscillations  as  a  result  of 
wave  dispersion  in  the  incident  bar  (figures  3  and  5).  This  dispersion  effect  will  be  absent  for 
direct  impact.  Nevertheless,  there  will  be  some  slight  deceleration  of  the  striker  bar  as  well  as 
(possibly  large)  oscillations  in  velocity  throughout  the  specimen  at  early  times  due  to  multiple 
reflections  from  the  specimen-bar  interfaces.  Unfortunately,  we  do  not  have  any  plots  of  the 
velocity  history  for  this  case. 
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In  the  SHPB  simulation,  gaps  formed  at  the  S-IB  and  S-TB  interfaces  about  25  ps  after  the 
arrival  of  the  loading  wave  at  the  S-IB  interface.  In  the  direct  impact  simulation,  gaps  formed 
much  sooner  at  both  the  specimen- striker  bar  and  S-TB  interfaces.  In  fact,  by  6  ps  after  impact, 
the  entire  face  of  the  specimen  had  separated  from  the  transmission  bar.  This  complete 
separation  continued  for  the  short  duration  of  the  simulation,  which  was  stopped  at  about  13  ps 
after  impact,  corresponding  to  a  compressive  axial  strain  of  about  3%.  Prior  to  gap  formation, 
the  peak  compressive  axial  stress  in  the  specimen  was  about  8  times  that  for  the  SHPB 
simulation.  Also,  tensile  axial  stresses  on  the  centerline  were  first  observed  in  the  transmission 
bar,  analogous  to  the  SHPB  simulation  (section  6.3). 

Our  goal  in  doing  these  direct  impact  simulations  was  not  to  study  the  direct  impact  compression 
test  per  se,  but  rather  to  examine  the  effect  of  a  different  loading  history  on  the  formation  of 
gaps.  However,  the  results  described  above  call  into  question  the  use  of  direct  impact  on  soft, 
nearly  incompressible  specimens.  Indeed,  whereas  inertial  effects  at  early  times  in  SHPB  tests 
can  be  reduced  (though  not  eliminated)  with  better  pulse  shaping,  for  a  direct  impact  test  there  is 
no  way  to  reduce  the  acceleration  impulse  imparted  to  the  specimen. 

C-l.l  Lateral  Constraints 

The  lateral  surfaces  of  the  specimen  and  the  bars  are  unconstrained  and  hence  stress-free  in  both 
the  SHPB  and  direct  impact  compression  tests.  We  performed  two  additional  direct  impact 
simulations  to  examine  the  effect  of  these  free  surfaces  on  gap  formation.  The  conditions  in  both 
simulations  were  as  described  above  for  the  direct  impact  test,  with  the  exceptions  noted  below. 

In  the  first  simulation,  the  lateral  surfaces  of  the  bars  and  the  specimen  were  constrained  to 
prevent  any  radial  motion,  whereas  the  axial  component  of  the  motion  at  these  surfaces  was 
unconstrained  as  before.  No  gaps  formed  in  this  case.  Axial  stress  contours  at  10  ps  revealed  a 
uniform  compressive  stress  state  throughout  the  specimen  and  the  bars  (in  the  neighborhood  of 
the  specimen),  with  the  exception  of  a  small  region  of  tensile  stress  adjacent  to  the  portion  of  the 
faces  of  the  bars  overhanging  the  specimen.  These  results  are  not  surprising.  Indeed,  under 
these  conditions,  the  specimen  essentially  undergoes  a  compressive,  uniaxial  strain,  just  as  in  a 
normal  plate  impact  test.  Since  the  radial  motion  throughout  the  specimen  is  eliminated  for  the 
most  part,  so  are  the  associated  radial  inertia  effects. 

In  the  second  simulation,  the  lateral  constraint  on  the  specimen  was  removed,  but  the  constraint 
on  the  bars  was  retained.  Gaps  formed  at  both  interfaces.  The  results  were  similar  to  those  for 
the  unconstrained  direct  impact  simulation  described  previously. 

These  two  simulations  support  the  conclusion  that  radial  inertia  effects  are  responsible  for  gap 
formation  in  the  direct  impact  tests  as  well. 
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C-2  Annular  Specimens 

All  of  the  results  presented  up  to  this  point  have  been  for  solid  (disc  shaped)  specimens.  Annular 
(washer  shaped)  specimens  have  also  been  used  in  SHPB  tests  on  soft  materials  in  an  effort  to 
reduce  inertial  effects  (5,  6).  The  motivation  for  hollowing  out  the  specimen  comes  from 
approximate  theoretical  estimates  for  the  inertially  induced  radial  stress  in  solid  specimens, 
which  predict  that  the  radial  stress  is  largest  on  the  centerline  and  decreases  parabolically  to  zero 
at  the  lateral  (stress-free)  surface  (6,  7,  26,  27,  32,  33).  Furthermore,  numerical  simulations 
reveal  that  in  situations  where  gaps  do  not  form  but  the  axial  acceleration  is  sufficiently  large, 
the  pressure  at  the  center  of  a  solid  specimen  may  undergo  large  oscillations  (32,  34).  Therefore, 
removing  the  center  portion  of  the  specimen  might  reduce  these  inertial  effects,  resulting  in  a 
stress  state  that  is  closer  to  uniaxial.  Experimental  results  indicate  that  this  works  for  early  times 
(5,  6)5 2 

We  performed  several  numerical  simulations  on  an  annular  specimen  with  the  same  initial 
thickness  ( Ls  =  1.45  mm)  and  initial  outer  radius  (Rs  =  6.35  mm)  as  the  solid  specimens 
discussed  previously.  The  inner  radius  of  the  annular  specimen  was  taken  to  be  2.638  mm.  As 
before,  we  considered  two  values  for  the  radius  of  the  pressure  bars,  Rb  =  12.8  mm  and 
Rb  =  9.5  mm.  The  dimensions  of  the  annular  specimen  and  the  smaller  bar  radius  were  chosen  to 
agree  with  those  used  in  the  experimental  study  on  ballistic  gelatin  by  Moy  et  al.  (5). 53  The 
compressible  Mooney-Rivlin  model  was  used  for  all  annular  specimen  simulations. 

For  the  larger  diameter  pressure  bars  and  a  25-ps  initial  rise  time  for  the  loading  wave,  the  total 
force  acting  on  a  cross-section  of  the  incident  and  transmission  bars  was  measured  at  an  axial 
location  of  2.52  mm  from  the  specimen-bar  interfaces  (cf.  section  7.3).  The  force  history  in  the 
transmission  bar  was  qualitatively  similar  to  that  for  the  solid  specimen  (figure  30),  although  the 
negative  forces  (indicating  gap  opening)  persisted  for  a  shorter  time.  The  force  history  in  the 
incident  bar  oscillated  between  positive  and  negative  from  157  ps  to  the  end  of  the  simulation  at 
340  ps  (43%  nominal  axial  strain),  indicating  repeated  gap  opening  and  closure. 

For  the  smaller  diameter  pressure  bars  and  a  25-ps  initial  rise  time,  no  gaps  were  observed 
except  momentarily  (at  about  7%  axial  strain)  at  the  S-IB  interface  over  a  small  region  midway 
between  the  inner  and  outer  surfaces. 

For  the  smaller  diameter  pressure  bars  and  a  1-ps  initial  rise  time  for  the  loading  wave,  a 
substantial  gap  was  observed  at  the  S-IB  interface  but  not  the  S-TB  interface.  At  t  =  182  ps  and 
an  axial  strain  of  7%,  this  gap  extends  over  %  of  the  S-IB  interface,  as  seen  in  figure  C-l.  This 
figure  also  reveals  bulging  of  the  annular  specimen  at  the  outer  (top)  and  inner  (bottom)  stress- 
free  lateral  surfaces,  near  the  S-IB  interface.  The  centerline  of  the  specimen  and  bars  lies  outside 

However,  theoretical  and  computational  results  demonstrate  that  the  presence  of  a  stress-free  inner  surface  leads  to 
additional  inertial  effects  at  later  times  (32,  36).  This  complication  will  not  be  discussed  further  here. 

53  Actually,  the  inner  radius  of  the  annular  specimen  used  here  is  1%  smaller  than  the  value  used  in  Moy  et  al.  (5). 
The  reason  for  this  slight  difference  is  that  it  allowed  us  to  use  the  same  mesh  size  for  the  solid  and  annular 
specimens. 
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of  (below)  the  figure.  Figure  C-l  also  shows  contours  of  the  pressure.  In  order  to  reveal  more 
clearly  the  regions  of  negative  (tensile)  pressure,  any  point  with  zero  or  positive  pressure 
(compression)  is  shaded  red;  any  color  other  than  red  indicates  a  tensile  state.  The  figure  reveals 
a  region  of  tensile  pressure  in  the  specimen  bordering  the  gap  at  the  S-IB  interface. 


Figure  C-l.  Pressure  contours  in  the  vicinity  of  an  annular  specimen  at  7%  axial  strain.  Any 

positive  pressure  (compression)  is  shaded  red;  other  colors  indicate  negative  (tensile) 
pressure  in  units  of  GPa.  The  loading  wave  had  a  1-ps  initial  rise  time.  Observe  the 
gap  along  most  of  the  S-IB  interface. 
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C-3  Simulations  with  the  PRESTO  Code 


With  the  exception  of  the  results  to  be  discussed  in  this  section,  all  of  the  numerical  simulations 
in  this  report  used  the  commercial  finite-element  code  LS-DYNA  (8,  9)  in  the  2-D  axisymmetric 
mode.  In  an  effort  to  verify  that  the  formation  of  gaps  at  the  specimen-bar  interfaces  is  not  some 
artifact  of  this  particular  code,  Bryan  Love  (37)  performed  analogous  numerical  simulations  with 
PRESTO,  a  finite-element  code  from  Sandia  National  Laboratories.  All  of  the  PRESTO 
simulations  were  fully  three-dimensional,  which  necessitated  larger  element  sizes  than  those 
used  in  the  LS-DYNA  simulations.  The  two  codes  also  have  different  contact  algorithms. 

A  compressible  version  of  the  Mooney-Rivlin  model  was  used  for  the  specimen,  with  the  same 
calibration  used  in  the  LS-DYNA  simulations  (see  section  3).54  Both  solid  and  annular 
specimens  were  considered;  these  had  the  same  dimensions  as  in  the  LS-DYNA  simulations. 

The  diameter  of  the  aluminum  pressure  bars  used  in  these  simulations  was  1  inch  (25.4  mm), 
which  is  essentially  the  larger  of  the  two  diameters  considered  previously.  The  pressure  bars 
were  much  longer  (2438  mm)  than  in  the  LS-DYNA  simulations,  corresponding  to  the  bar 
lengths  in  an  actual  Hopkinson  bar  setup  (25).  The  loading  condition  in  the  PRESTO 
simulations  involved  an  aluminum  striker  bar  (of  the  same  diameter)  impacting  the  far  end  of  the 
incident  bar,  as  in  a  real  SHPB  test.  However,  no  pulse  shaper  was  used  between  the  striker  and 
incident  bars.  A  few  direct  impact  simulations  were  also  performed. 

The  results  of  these  PRESTO  simulations  confirmed  the  findings  of  the  previous  LS-DYNA 
simulations  with  regard  to  the  formation  of  gaps  at  the  specimen-bar  interfaces.  The  fact  that  the 
PRESTO  simulations  were  fully  three-dimensional  rules  out  the  possibility  that  gap  formation  in 
the  LS-DYNA  simulations  was  in  some  way  associated  with  issues  at  the  centerline  in  the  2-D 
axisymmetric  mode. 


However,  the  PRESTO  version  of  the  compressible  Mooney-Rivlin  model  does  not  suffer  from  the  unrealistic  behavior  of 
the  pressure  at  very  large  volumetric  compression  and  high  Poisson’s  ratio  that  is  exhibited  by  the  LS-DYNA  version 
(see  appendix  A-3). 
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List  of  Symbols,  Abbreviations,  and  Acronyms 


Ai,  A 2 
Cl,  Cg,  Ce 

E 

ez 

Fz 

G 

IB 

J 

K 

L 

Em 

Ls 

Ljb 

P 

r 

R 

Rb 

Rs 

S 

SFACT 

SHPB 

S-IB 

S-TB 

t 

TB 

tR 


constants  in  the  Mooney-Rivlin  model 
longitudinal,  shear,  and  bar  wave  speeds 
Young’s  modulus 

nominal  axial  strain  in  the  specimen 

total  force  on  a  cross-section  of  the  transmission  bar 

shear  modulus 

incident  bar 

Jacobian  of  the  deformation 

bulk  modulus 

longitudinal  modulus 

length  of  the  incident  bar 

initial  length  (thickness)  of  the  specimen 

length  of  the  transmission  bar 

pressure 

radial  coordinate  of  the  deformed  material 
radial  coordinate  of  the  undeformed  material 
radius  of  the  incident  and  transmission  bars 
initial  radius  of  the  specimen 
specimen 

scale  factor  for  the  penalty  force  stiffness  in  the  contact  algorithm 

split  Hopkinson  pressure  bar 

specimen-incident  bar 

specimen-transmission  bar 

time 

transmission  bar 

rise  time  of  the  axial  velocity  prescribed  at  the  end  of  the  incident  bar 


104 


Wib,  «s-ib 

WTB,  Ms-TB 

Vo 

VDC 

Vib,  Vtb 
Vz 
z 
Z 

Sib,  Ss-ib 
STB,  Ss-TB 


As-ib 

As-tb 

V 

^-zi  ^ '/>  M, 

G 

p 

G zz ,  ®rn  ^00 


axial  displacements  of  the  incident  bar  and  specimen  at  the  S-IB  interface 
axial  displacements  of  the  transmission  bar  and  specimen  at  the  S-TB  interface 
plateau  value  of  the  axial  velocity  prescribed  at  the  end  of  the  incident  bar 
viscous  damping  coefficient  in  the  contact  algorithm 

axial  velocity  of  the  incident  and  transmission  bars  at  the  specimen  interface 

axial  velocity  of  the  incident  bar 

axial  coordinate  of  the  deformed  material 

axial  coordinate  of  the  undeformed  material 

deformed  axial  coordinates  of  the  incident  bar  and  specimen  at  the  S-IB  interface 

deformed  axial  coordinates  of  the  transmission  bar  and  specimen  at  the  S-TB 
interface 

gap  size  at  the  specimen-incident  bar  interface 
gap  size  at  the  specimen-transmission  bar  interface 
Poisson’s  ratio 

deformed  length  (thickness)  of  the  specimen 
axial,  radial,  and  hoop  stretches  in  the  specimen 
Cauchy  (true)  stress  tensor 
initial  density 

axial ,  radial,  and  hoop  components  of  the  Cauchy  stress 
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